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SUMMARY 


Results of a study towards the development of flutter modules applicahle 
to automated structural design of advanced aircraft configurations, such as a 
supersonic transport, are presented. In this study automated structiiral 
design is restricted to automated sizing of the elements of a given structural 
model. It includes a flutter optimization procedure; i.e., a procedure for 
arriving at a structure with minimum mass for satisfying flutter constraints. 
Methods of solving the flutter equation and computing the generalized aero- 
dynamic force coefficients in the repetitive analysis environment of a flutter 
optimization procedure have been studied and recommended approaches are pre- 
sented. Five approaches to flutter optimization are explained in detail and 
compared. An approach to flutter optimization incorporating some of the 
methods discussed is presented. Problems related to flutter optimization in 
a realistic design environment are discussed and an Integrated approach to 
the entire flutter task is presented. Recommendations for further investiga- 
tions are made. Results of numerical evaluations, applying the five methods 
of flutter optimization to the same design task, are presented. 
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STUDY OF FLUTTER RELATED COMPUTATIONAL 


PROCEDURES FOR MINIMUM WEIGHT STRUCTURAL 
SIZING OF ADVANCED AIRCRAFT 
R. F. O'Connell, H. J. Hassig and N. A. Radovcich 
Lockheed-California Company 
Burbank, California 


1 . INTRODUCTION 
1 . 1 General 


One of the factors contributing to the profitability of an airplane is its 
payload/range capability. Given certain safety and performance requirements, 
there is a direct trade-off between structural weight and payload, and it is 
the ideal of each airplane designer to reduce structural weight. Although the 
ideal minimum weight design may be expensive to produce, overshadowing any 
payload/range gains , it provides a good starting point for a practicable 
design and a good basis for comparing different designs. 

Structural weight minimization, of course, is not a new idea. It is one 
of the airplane designer's most critical tasks. It now has come to the fore- 
front as a result of two developments. 

First, it has become evident to the structural design engineer that the 
combination of finite element modeling, high speed computer capacity, and 
mathematical techniques maJtes it practicable to do detailed structural 
synthesis aimed at minimizing weight. 

Second, the need for a comprehensive and detailed approach to structural 
design optimization has significantly increased with the advent of the super- 
sonic transport . This follows from the fact that for a supersonic transport 
the return in terms of payload/range per pound of structural weight saved is 
considerably larger than for a subsonic transport. For instance, a one per- 
cent structural weight saving on a typical subsonic transport might result in 
an increased payload capability of one to two percent; on an arrow wing super- 
sonic transport , recently studied by the Lockheed-California Company , a one 
percent structural weight saving resulted in a four percent increase in pay- 
load capability for the design range. 

The subject of this report is flutter optimization; i.e., structural 
weight minimization with flutter constraints. The need for a systematic, 
possibly automated, approach to flutter optimization also has increased sig- 
nificantly. Subsonic transports, as they are known, and transonic transports. 



as shown in artist's sketches, can he represented by simple, beam-type 
structural models that are satisfactory for optimization with flutter con- 
straints . Flutter optimization for such designs can be done, and has been 
done, with available methods. The supersonic transports that are flying and 
those being studied, however, all have lifting surfaces that cannot be repre- 
sented satisfactorily by simple beam- type structural models. This fact alone 
makes the task of flutter optimization an order of magnitude more complicated. 

Although ad hoc approaches to flutter optimization still could lead to a 
satisfactorily optimized supersonic design, refined methods that take full 
advantage of the capabilities of the present computers, in regard to automation 
as well as interaction with the engineer, become attractive and possibly 
mandatory. This is especially true in view of the rapidly increasing capa- 
bility for fast analysis and synthesis in the areas of structural modeling and 
analysis, stress optimization, and performance analysis supported by improved 
configuration control. Flutter optimization must keep abreast of these devel- 
opments. A balanced improvement in capability in all disciplines will make 
possible, within a practicable time span, true in-depth comparisons between a 
large number of candidate designs. 

The preceding paragraphs present generally well known justification for 
a concerted effort in improving methods of structural optimization with flut- 
ter constraints. Work performed dixring the subject study is part of such an 
effort . 

Work towards the goal of a generally available automated or semi -automated 
structural optimization system, that includes items such as optimization for 
stress, flutter and controllability, multiple flutter speed and modal damping 
constraints, is still in a state of development. The present study has con- 
tributed to this goal in the following areas. Methods of computing the aero- 
dynamics parameters to be used in a flutter optimization program have been 
compared in detail with respect to characteristics which are independent of a 
specific aerodynamics theory. A method for efficiently and reliably solving 
the flutter equation for roots of interest in a flutter optimization module 
has been developed. Five methods of flutter optimization have been compared 
in detail and the mechanics of the optimization process have been examined; 
nimierlcal examples with all five methods have been generated for the same 
aircraft design. Recommendations for further study and for the design of a 
flutter optimization module have been made. 

The principal results of this study are presented in this report. Back- 
ground discussions and supporting material are presented in a companion report 
(Reference l). 
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1.2 Objectives of Study 


The objectives of this study are: 

1. To survey and evaluate methods of representing unsteady aerodynamics 
parameters and make recommendations for a general, accurate and efficient 
formulation that minimizes the computational effort during the optimiza- 
tion process. The assessment of aerodynamics theories, however, falls 
outside the scope of this study. 

2. To survey and evaluate methods of determining the flutter characteristics 
and make recommendations for a method that is reliable and efficient, and 
suitable for the optimization process. 

3. To evaluate and compare a number of methods of structural optimization 
with flutter constraints and make recommendations for further evaluation 
in a realistic design environment. 

4. To make preliminary recommendations for the design of a flutter 
optimization module. 


2. OVERVIEW OF THE FLUTTER OPTIMIZATION TASK 


Structural optimization with flutter constraints is both an extension of 
the structural optimization task related to strength and an extension of the 
flutter analysis task. Being an extension of two tasks that traditionally are 
considered to belong to different disciplines, flutter optimization must take 
into account requirements of both disciplines. Structural optimization 
requires that a structural model is used that incorporates sufficient struc- 
tural detail, in terms of distribution of structural material, to aid the 
designer in defining hardware. Similarly the flutter analysis that is incor- 
porated in the optimization process must be of an accuracy comparable to that 
used outside of flutter optimization. The latter refers to methods of 
representing the unsteady aerodynamics and methods of solving the flutter 
equation, since the more detailed a structural model is, the more accurate, 
from an idealized theoretical point of view, is the flutter analysis. From a 
practical point of view, structural sizing for strength requires more detail 
in the structural model than is required for adequate prediction of flutter 
characteristics . 

Thus, the flutter optimization task starts with the definition of the 
structural model. This is one of the most crucial aspects of flutter optimiza- 
tion, and it involves a serious conflict between simplicity of approach and 
computer cost. Present computer technology, or methods of structural analysis, 
or both, may not permit a structiiral model with sufficient detail for a stress 
analysis to be used in flutter optimization; computer cost could be exorbitant 
due to the repetition of operations during the design process. Section T-1 
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deals with this problem in more detail. Suffice here that associated with the 
choice of structural model is the selection of a practicable mmiber of degrees 
of freedom for the vibration analysis that has to provide the modes for the 
modal reduction of the flutter equation, which is usually required to limit 
computer cost. If the degrees of freedom for the vibration analysis are a 
subset of the degrees of freedom of the structural model, complications arise 
if a nonlinear relationship between the stiffness matrix and the design 
variables results (see Section T-l)- 

Without serious restriction on scope or accuracy of the analysis the mass 
matrix can be assumed to be a linear function of design variables. It is the 
sum of a basic matrix and as many elementary matrices as there are design 
variables associated with a mass change, each proportional to a design 
variable . 

During the flutter optimization there is repeated need for determining 
roots of the flutter equation, each time that a structure that has undergone 
a resizing since the previous solution of the flutter equation. For many, if 
not all, of these solutions a remodalization is necessary based on vibration 
modes of the current configuration. It is found that for the optimization 
process to provide reliable, converged results, consistent with the capability 
of the structural model, more modal degrees of freedom are required in the 
flutter equation than for a routine flutter analysis (see Section i+). 

Incorporation of state-of-the-art level aerodynamics in the flutter 
optimization process does not provide significant problems beyond those 
encountered in the usual flutter analysis. For a given external geometry the 
basic aerodynamics formulation is invariant with structural changes. The 
repetitive formation of generalized aerodynamic forces for successive, updated 
modalization of the flutter equation is simple and relatively inexpensive. 

In view of the objectives and the scope of the present program, this 
report devotes major sections to important aspects of the flutter optimization 
procedure. Section 3 deals with the solution of the flutter equation. Sec- 
tion h deals with modalization. Section 5 presents part of considerable work 
devoted to the aerodynamics, with the remainder being presented in Reference 1. 
In Section 6, methods of optimization for flutter, evaluated during this study, 
are discussed. Numerical results obtained by applying these methods to a 
simplified optimization task are presented in Appendix A. 

Against the background provided by these sections. Section 7 presents 
discussions of several additional problems and considerations that need to be 
studied in order to choose a rational approach for formulating a flutter 
optimization module. 

In Section 8, computational aspects of the complete flutter task are 
delineated. This task includes flutter analysis as well as structural syn- 
thesis of a design that satisfies the flutter requirements. 

Section 9 summarizes the conclusions of the present study and presents 
recommendations for future work. 
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3. SOLUTION OF THE FLUTTER EQUATION 
3.1 The Generalized Flutter Equation 


When using the k method the flutter equation can be written as : 





0 


(3.1) 


One of several possible methods of solving this equation is to determine the 

characteristic value X. = for several values of the reduced frequency 

UiC f \ 

k=— ^ , keeping all other quantities in the equation constant (Reference 2) . 

In the p-k method the flutter equation is 



(3.2) 


and solutions p=(y+i)k are sought for selected combinations of values of 
V and p (Reference 3). The p-k method formulation is convenient for the 
inclusion of viscous damping and control system transfer functions. This is 
accomplished by -writing; 


i- p2 p + (l+ig) [k] --|-pV^[A(ik)] -SH.(p)[d.] {q} = 


(3.3) 


where H.(p), j=l,2.., represents transfer functions of the control system 
(D 


and 




relates the control system displacements to the structural dis- 


placements ; [D] is a viscous damping matrix (Reference 3) 


A further generalization of the flutter equation can be made by making 
the stiffness matrix and the inertia matrix functions of design variables 


m^, which is the standard proced-ure for structural optimization. In addi- 


tion, other quantities, such as 


H • W 


as well as transfer function 


coefficients in H (p), may be made functions of design variables. 
J 
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Equation (3.3) implies that the determinant of the square matrix on the 
left hand side is zero and thus, in a very general form, the characteristic 
equation corresponding to the flutter equation can be written as: 

(3.U) 


D is called the flutter determinant . For arbitrary values of the variables 
it has a complex value. Thus equation (3.4) represents two equations and, in 
principle, can be solved for two unknowns for given values of the other 
variables . 

Letting 7=0 and solving for g and V corresponds to the traditional 
k method of solving the flutter equation. Solving for 7 and k corres- 
ponds to the p-k method. Letting 7=0 and solving for k and V leads 
directly to the flutter speed for a given value of the structural damping, g. 

Solving equation (3.4) for k and one of the design variables, assuming 
all other variables fixed, is a new use of the flutter equation. It is called 
Incremented Flutter Analysis (References 4 and l), applications of which are 
included in Sections 6.2.3 and 6.6. 


d| (7+i)k,g,V,p,m^| = 0 


3.2 Types of Solution Sought 


A flutter analysis in the traditional sense is the determination of the 
flutter characteristics of a given structure. It includes the calculation 
of any flutter speed that may occur at speeds up to or somewhat beyond a 
speed corresponding to the required flutter margin. It also includes the 
gaining of insight in the variation of frequency and aerodynamic damping at 
speeds below the flutter speed for several in-flight vibration modes of 
interest. Consequently, sufficient modal solutions are obtained for the con- 
struction of f-g-V diagrams (Figure 3-l) for several flight conditions. 


A procedure for structural optimization with flutter constraints will 
most likely start with such a survey-type analysis. However, during the 
process of repeated resizing, leading to the optimum design, there is no need 
for determining complete f-g-V diagrams at each resizing; only point solu- 
tions are required. Point solutions found in the literature are of two types: 
l) directly solving for the flutter speed (the combination k,V in equa- 
tion (3.4)), and 2) determining the value of one design variable necessaiy to 


satisfy a given flutter speed constraint (the combination 


k, m . 


tion (3.4) where 


m. 

J 


is one of the design variables 



in equa- 
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One additional type of point solution has been formulated during this 
contract , resulting in the determination of the minimum damping point of an 
in-flight mode. Such a point, if it exists, is of interest if the minimiam 
damping point lies within the speed range considered. The associated mode is 
called a hump mode (See Figure 3-1 ). This point solution requires the solution 
of equation (3.^) and the equation 



( 3 . 5 ) 


for the three unknowns k, 7 and V. Details of the formulation are given in 
Section 3.^. No numerical evaluation of the method has been made thus far. 

The following section deals mainly with methods of obtaining point 
solutions for the flutter speed. It should be kept in mind that when such 
solutions are needed in an optimization program a solution for a similar 
structural configuration is usually available as a first approximation to the 
required solution. 


3.3 Methods of Obtaining Point Solutions 


Several methods for obtaining point solutions have been considered and 
evaluated to various depths. Their apparent efficiency, in terms of computa- 
tional effort, is an important part of the evaluation. However, the degree 
of certainty with which a desired solution can be found is even more 
important . 

The latter consideration refers to convergence problems and to problems 
associated with relating modal solutions at one value of V or k to modal 
solutions at another value of V or k. 

The results of evaluations of the following methods are presented: 

• Bhatia method (Reference 5) 

• Phoa - Boeing method (Reference 6) 

• Lockheed's Program 165 (p-k method. Reference 3) 

• Desmarais-Bennett method (Reference 7) 

• Two Dimensional Regula Falsi and Newton Raphson (Reference 8) 

3 . 3.1 Bhatia Method - In Reference 5 Bhatia presents a method of solving 
directly for the flutter speed. Nimierical evaluations of the method have 
been performed using data from the arrow wing study that Lockheed has con- 
ducted under contract NASl-12288. 

In Bhatia' s method, which is based on the k-method approach, the struc- 
tural damping, g, required for neutral stability is computed as a function 

of 1/k = — — by means of a Laguerre type extrapolation. It is an iterative 
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method that is initiated hy choosing a trial value l/k and computing the 
associated value of g and its first and second derivative with respect to 
l/k. The Laguerre extrapolation leads to a first approximation of the value 
of l/k for which g=0... The process is repeated for this new value of l/k 
until convergence is reached. 

The method as presently programmed uses only aerodynamic matrices at 
preselected values of k, requiring a large number of preselected k values. 
The method also requires inputting the first and second derivatives of all 
aerodynamics matrices with respect to l/k. The method could be improved by 
using interpolation with respect to k to determine the aerodynamics matrix 
and its derivatives at arbitrary values of k from matrices given at a 
moderate number of preselected k values. Care must be taken that the 
interpolated results are defined -uniquely over the range of k of interest 
for a particular solution to prevent the solution from oscillating between 
two values ("hunting"). 

Nijmerical evaluations were performed as part of this study. Difficul- 
ties were encountered in tracking the proper mode and in converging on the 
lower flutter speed of a hump mode. There is uncertainty whether the program 
can be modified such that the proper modal solution is always obtained. 

At each step in the iteration towards the solution a characteristic 
value problem must be solved. This may prove to be costly in terms of CPU 
time . 


3.3.2 Phoa Method - In Reference 6, Phoa presents a formulation of the 
flutter equation from a controls theory point of view. Although it is 
recognized that controls theory could prove to be of assistance in inter- 
preting the flutter phenomenon, in the case of Reference 6 it leads to an 
equation that is essentially the same as equation (3.4). Phoa's method, based 
on the k-method approach, is in use at the Boeing Company. Discussions with 
Boeing personnel indicate that in the actual application the equation 
{d(o),V). -l} = D(oj,V) = -1 is solved for V and w. 

The solution is accomplished in two steps. Constant velocity lines in 
the complex plane representing D(co,V) are intersected with the real axis. 

The values of the real parts at the intersections, as a function of the 
velocity, are used to determine an estimate of the velocity for which the 
real part of D(w,V) equals -^1. In an iterative process the accuracy of 
the solution is improved. 

In numerical evaluation of this approach it was shown that the constant 
velocity lines may have two intersections -with the real D(oj,V) axis; this 
can be a source of problems (Figure 3-2). 

The method is a sequence of two interpolations requiring many determinant 
evaluations. It is expected that very few, possibly not more than two or 
three, steps in the Iteration process are required. 
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3.3.3 Lockheed Program 165 - Program I 65 of Lockheed's Flutter and Matrix 
Algebra System (FAMAS) is based on the p-k method approach. It is designed 
to generate many point solutions, associated with in-flight modes, so that com- 
plete f-g-V diagrams can be constructed. The program solves equation C3^3) 
for p=(Y+i)k given an initial trial solution. For determining the flutter 
speed, y is evaluated at several values of the velocity. Flutter occurs at 
the speed for which Y=0 . 

The program has been used successfully in nonautomated numerical evalua- 
tions during this contract. Automation should be relatively simple and cotild 
be based on the following steps. At the estimated flutter speed an estimated 
frequency is used to start the process. Both are obtained from the solution 
for a previous structiiral configiaration. Determinant iteration (see 
Reference 3) will lead to the actual value of Y at the estimated flutter 
speed. At a slightly perturbed velocity, using the damping and frequency 
already found as trials, Y is again evaluated. The two pairs of V and Y 
values thus found are used to initiate a One-Dimensional Regula Falsi pro- 
cedure that leads to a value of V for which Y=0. The approach is expected 
to be quite efficient, except for the problem of assuring that subsequent 
solutions belong to the same in-flight mode. 

3 . 3 .U Desmarais-Bennett Method - Reference T presents a fast and economical 
automated procediire to generate f-g-V diagrams, including the proper con- 
nection of point solutions of the flutter equation. The procedure is based 
on the k-method approach . 

Reference 7 shows that the method is quite powerful in properly con- 
necting point solutions. The sample cases in Reference 7» however, are 
obtained by partial deflation of the flutter determinant after each modal 
solution is found. Thus using this method would require solving for more 
roots than are of interest if only the flutter speed is required. Or, 
alternatively, if only the root of interest is determined, there is uncer- 
tainty whether the method will be as successful in following modal solutions 
as shown in Reference 7. 

Application of this method to directly solving for the flutter speed 
coTild be programmed according to the following procedure. 

The known solution for a base configuration is considered a reasonable 
estimate of the solution for a slightly modified configuration. Two k 
values, closely spaced according to the Desmarais-Bennett approach, are 
chosen such that flutter is expected to occur at a lower k value. Modal 
solutions at these two k values are obtained. The repeated sequence of 
linear extrapolation to the next k value and the Laguerre iteration 
described in Reference 7 is performed for the mode that is expected to give 
the flutter crossing and one or more additional modes on each side of this 
mode in the frequency spectrum. The additional modes are included to assiire 
that a flutter crossing is obtained, in the event that an error in Judgment 
is made in selecting the prime candidate mode for a flutter crossing. 
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The preceding conceptual evaluation defines the problems that need to be 
resolved when adjusting the Desmarais-Bennett method for use in a flutter 
optimization program and no numerical evaluation was considered necessary. 

3 . 3.5 Two-Dimensional Regula Falsi - The concept of solving the two equations 
implied by equation (3.4) for two unknowns is not new. However, using this 
concept for directly solving the flutter equation for the flutter speed is 
relatively new. The need for such a solution arose with the advent of struc- 
tural optimization with flutter speed constraints and, to the knowledge of 
the present authors, the first published record of solving directly for the 
flutter speed is Reference 9- 

In that Reference the Newton-Raphson approach is used in two dimensions 
to determine flutter speed and, as a byproduct, flutter frequency. The 
Newton-Raphson approach is based on determining the value of a function and 
its derivatives for an initial set of trial values and extrapolating linearly 
to an estimate of the solution. In Reference 9» the derivatives are deter- 
mined by a finite difference technique. The Two-Dimensional Regula Falsi 
approach uses three trial sets of the unknowns to construct two planes. The 
common point between those planes and the plane D(w,V)=0 defines the next 
estimate of the solution. 

Table 3-1 compares the essential characteristic of the two methods. In 
the Newton-Raphson method with analytical evaluation of the derivatives , the 
formation of two derivative matrices is time consming. In all methods the 
determinant evaluations are the most time consuming. Other operations, 
related to solving two linear equations with two unknowns, are trivial. Pro- 
visions to assure convergence are comparable for the two methods. Niomerical 
experience with the Two-Dimensional Regula Falsi has indicated that problems 
with convergence on a solution are more easily solved than with the Newton- 
Raphson approach. It is concluded that the Two-Dimensional Regula Falsi 
approach is the more preferable one of the two. 

It should be noted that both methods can be used for combinations of 
unknowns other than frequency and flutter speed. The Two-Dimensional Regula 
Falsi has been used successfully for solving for the value of one design 
variable, required to meet a given flutter speed, and the associated frequency. 
The method does not require the computation of derivatives . No interpolation 
or extrapolation of converged solutions is required, unless nonconvergence is 
encountered and an intermediate configuration is analyzed to assist in obtain- 
ing a better initial estimate of the solution for the configuration for which 
the original nonconvergence occurred. Finally, the solution sought is a 
combination of real values of the unknowns, rather than a series of complex 
modal solutions associated with in-flight modes. The equivalent of converging 
on the wrong mode, as may occur in seeking modal solutions, usually leads to 
nonconvergence, and a recovery procedure that is described in Reference 1. 

Thus mode switching to a non-flutter mode does not occur or, at worst , leads 
to nonconvergence. The chance of converging on the wrong flutter speed and 
frequency would seem to be quite small in view of the relatively small number 
of solutions within the region of interest of the unknowns. It has never 
occurred in the many test cases that have been run during this study. 
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TABLE 3-1. COMPARISON OF NEWTON-RAPHSON METHOD AND 
TWO-DIMENSIONAL REGULA FALSI METHOD 


Operation 

Newton-Raphson 

Two- 

Dimensional 
Regula Falsi 

Analytical 

Derivatives 

Finite 

Differences 

Number of initial estimates 

1 

1 

3 

Interpolation of aerodynamics 
matrix required? 

Yes 

Yes 

Yes 

Derivative of aerodynamics 
matrix required? 

Yes 

No 

No 

Formation of derivative matrices 
required? 

Yes 

Trivial 

No 

Number of complex 
determinant 
evaluations per 
iterative step 

First 

step 

Value of 
determinant 

1 

3 

3 

Derivative 

2 

0 

0 

Each 

fol- 

lowing 

step 

Value of 
determinant 

1 

3 

1 

Derivative 

2 


0 


3-3.6 Conclusion - On the basis of overall engineering evaluation, supported 
by ntunerical experience vith all methods except the Desmarais-Bennett and 
Newton-Raphson approach, the Two-Dimensional Regula Falsi approach was con- 
sidered most promising and chosen for further development (see Reference l) . 


3.^ Minimum Damping in Hump Mode 


Sufficient modal damping within the speed envelope can be assured by 
requiring sufficient damping in all modes at "all" speeds below the minimum 
required flutter speed or by requiring that the minimum damping in each mode, 
in so far as it occurs below the minimimi required flutter speed, is equal to 
or larger than a given value . 

To initiate exploration of the latter approach a method to determine the 
minimum damping in a hump mode was formulated. 
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The point of minim-um damping in the hump mode is defined by the condition 


= 0, where Y defines the real part of the flutter root, p=(7+i)k, in 
d V 

terms of the reduced frequency k. The quantity 7 is a form of the 
logarithmic increment: 

(3.6) 

n 

where a and a , are amplitudes of successive cycles . 
n n+1 

An expression for is found as follows . 

a V 

Consider the p-k method formulation of the flutter equation (equa- 
tion (3.2)) and take the derivative with respect to V: 


21 [„] , 2 vi [„] p |£ - pv [a( I k)] - i pv" [A{lk)] -fij j , } 




(3.T) 


37 


Choose a velocity for which is estimated to be equal to zero. 

1 <7v 

The solution of the flutter equation at is: p = p^, jqj = jq^j and the 

characteristic vector of the transposed equation: jrj = • Substituting 

this solution into equation (3-7) gives: 


[«]{%} 


2 ^1 

M M ^1 1 Pi + 2 — 

C 




-pVi[--J[A(ik)]l,p)tipv/ [rj [A'Ukp)] (3.8) 


where j^A'(ik^) J ~ |^A(ik)J evaluated at k = k^. 
With p = (7+i)k: 



Substituting equation (3*9) into equation (3.8) leads to a complex equa- 


tion, and thus two equations in the two unknowns and ^ from which 

d V dV 


(w), 


can be determined. 

The process can be repeated for V , leading to (■f^) . A one- 

d \dv/2 

dimensional Regula Falsi approach will lead to the value of V for which 



In the above approach two characteristic value problems must be solved 

97 

for detennining the first iterated value of V for = 0. Each following 

9 V 

step requires solution of one characteristic value problem. 


It should be noted that damping versus speed curves may be rather flat 
and for practical purposes a converged value of V may not define a converged 
value of V. This causes no problem since the most likely application of 
these procedures is in connection with an inequality constraint such as: 


hxamp top “ max allowed 


( 3 . 10 ) 


Determining the minimum damping in a hump mode can be combined with 
solving for the value of a design variable satisfying the constraint: 


7 =7 =7 (3.11) 

hump top max allowed 

For V = and 7=7, equation (3.7) is solved for k and the value 

of the design variable m^ . Then equations (3.9) and (3.8) are used to compute 

^ as before. In general ~ ^ 0. and a one-dimensional Regula Falsi process 
dV dV 

is initiated by repeating the process for another chosen value V = V^- 

Numerical evaluation of the approaches outlined could not be accomplished 
within the scope of this study. 
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3 . 5 Recommendation 


The two-dimensional Regula Falsi procedure is recommended for inclusion 
in the Flutter Optimization Module for providing point solutions of the 
flutter equation. . The procedure is more direct than any of the other pro- 
cedures considered. It aims at roots of the flutter equations, either flutter 
speed and frequency or design variable and frequency, of which for every 
flight condition there are considerably fewer present than there are in-flight 
modes represented in the problem formulation. As a result, convergence on the 
wrong root would seem to be less likely than when modal solutions are sought. 
That the same procedure can be used for solving for different pairs of 
unknowns is considered an added advantage. In addition, it is equally appli- 
cable to the p-, the k- and the p-k method of formulating the flutter 
equation. A preliminary program is available that has shown good convergence 
behavior under a wide variety of input data. 


Since it seems likely that the capability of directly solving for the 


point for which 



will be a factor in developing methods of flutter 


optimization, numerical test cases should be conducted to evaluate the methods 
related to determining the minimum damping in the hump mode. The results may 
influence the development of methods of optimization that take into account 
damping constraints. 


1+ . MODALIZATION 
4 . 1 General 


Modalization is the reduction of the number of degrees of freedom by 
establishing modes of displacement in which the original degrees of freedom 
(usually point displacements) have a fixed relation to each other. 

Let define a relation between the discrete structural displace- 

ments z. The arbitrary column matrix of displacements {z} can then be 
approximated , by linear combination of several linearly independent columns 
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or in short notation: 


{z} = [z]{q} ( 1 ^. 2 ) 


The modalized flutter equation is: 

+ (l+ig)[K] - l/2pV^[A(ik)]j[^zj|q| = 0 (U.3) 

' Modalization is desirable whenever the total number of Initial degrees of 
freedom is sc large that solving the unmodalized equation becomes uneconomical, 
and is necessary if the number of initial degrees of freedom exceeds the 
capacity of the available computer program to solve the original characteris- 
tic value problem. Since, in general, the flutter equation is solved more 
frequently than the vibration equation and, in addition, the flutter equation 
must admit complex numbers, modalization is usually associated with the 
flutter equation. However, when using all the structural displacements of a 
detailed finite element structural model as degrees of freedom, modalization 
may be desirable or necessary for the vibration analysis as well. 

In any discussion of modalization, the type of modes and the number of 
modes to be used must be considered. When used in an optimization procedure, 
the question of "updating" must be considered. Updating in this context 
means the adjustment of the modes after resizing the structural elements in 
the course of the optimization procedure. These three aspects of modalization 
will be discussed separately in the following sections. 


U . 2 Types of Modes 


Before the advent of the high-speed computer, modalization (e.g., 
Rayleigh-Ritz method) was required even for vibration analyses. Relatively 
few and simple modes were used. With the increasing capacity of computers, 
the need for modalizing the vibration equation has all but disappeared. Thus, ■ 
present practice is to determine natural vibration modes of the entire airplane 
from an unmodalized vibration equation and to use a certain number of modes, 
associated with the lower range of natural frequencies, to reduce the order of 
the flutter equation. For special investigations, such as the inclusion of 
actual control-surface-actuator impedances, or the entire automatic control 
system, additional control surface modes may be necessary. 

In several instances in the literature Ce.g. , Reference lO), the use of 
component modes has been described. Component modes define the relations 
between discrete displacements of airplane components such as the wing or 
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fuselage, and are obtained by a vibration analysis in which only displacements 
of a particular component are used as degrees of freedom. Complications arise 
when the connections between components involve many structioral displacements. 
Reference 11 shows, with a simple beam as an example, that the unjudicious 
use of component modes can give inaccurate results for even the lowest fre- 
quency of the entire body. The use of component modes is only recommended 
for the determination of natural vibration modes of the complete vehicle, 
and then only if it is necessary to reduce the order of the vibration 
equations . 

Analytical modes, such as defined by polynomials and modes corresponding 
to static deflections, would obviate the need for repetitive vibration analyses 
during the optimization process if they are used as fixed modes. However, 
usually a considerably larger number of such modes is required, for the same 
accuracy of the flutter solution, than when vibration modes are used. No 
advantages off-setting that disadvantage have been encountered. 

Specifically, analytical modes have been suggested for efficient genera- 
tion of generalized aerodynamic force coefficients, as discussed in Section 5. 
The use of analytical modes may permit the analytical integration of the 
product of deflection and pressure modes. It makes it possible to compute 
invariant generalized aerodynamic force coefficients that can be combined 
linearly to form generalized aerodynamic forces for any arbitrary mode. To 
take advantage of this feature, however, the number of analytical modes must 
be large, since it must be adequate for a large number of stiffness and 
inertia distributions . The analytical modes thus can serve as reference modes 
that are the degrees of freedom for all vibration analyses from which a 
smaller number of vibration modes is obtained for use in flutter calculations. 
However, a large number of vibration modes of a basic configuration also can 
be used as reference modes and one would expect that fewer reference modes 
are needed if they are vibration modes than if they are analytical modes. 

It was thought that using the (complex) flutter mode of a base configura- 
tion might reduce the number of modes required for an adequate flutter solution 
of a modified configuration. Some preliminary work during this study was 
done, but was not carried far enough for any conclusion to be drawn. 


4.3 Number of Modes 

The number of modes used in the flutter equation is of importance for the 
accuracy of the computed flutter speed and flutter speed derivatives with 
respect to design variables. At present there seems to be no readily available 
general criterion for determining the member of modes needed for a desired 
accuracy. 

When trying to economize by restricting the number of modes to be used in 
flutter calculations, there is a need to frequently check whether the number of 


18 



modes is sufficient for acciirate prediction of the flutter speed for arbitrary- 
configurations. Thus there is an advantage in using flutter analysis pro- 
cedures that allow a large number of modes even if that raises the cost of 
each individual flutter solution. 

It has been pointed out in Reference 12, and it was confirmed by limited 
numerical analysis during this study, that more modes are needed for accurately- 
computing flutter speed derivatives than for computing flutter speeds . 

In deciding on the n-umber of modes the computer environment may be an 
Important factor to be Judged by the analyst in addition to the accuracy 
required. Even the method of computer cost appropriation may influence the 
decision. 


Updating of Modes 


As resizing steps accumulate during the optimization process, the vibra- 
tion modes of the initial configuration become less suited to accurately 
represent the revised structure. Ideally, therefore, after each resizing step 
a new vibration analysis should determine new modes for modalizing the flutter 
equation. The need for such updating is closely related to the number of 
modes used and the type and magnitude of structural changes incurred by the 
resizing. The use of a large number of modes tends to reduce the need for 
frequent updating. However, insufficient updating can cause the resizing 
steps to follow a zig-zag path that, in the extreme, may not converge. 

The physical explanation for this is the following. Let the optimization 
procedure indicate a local stiffening as the optimum resizing for resizing 
step j. Then the vibration modes for step j+1 would show a decrease in local 
deformation. If the vibration modes for step j are used for step j+1, the 
excess local deformation tends to reinforce and overestimate the beneficial 
effect of that local stiffening. Thus, in the absence of modal updating, 
material tends to be added where the first resizing step, with the modes used, 
indicates where it is most beneficial. 

Modal updating must not be confused with making the modal matrix a func- 
tion of the design variables. This aspect of modalizing was recently intro- 
duced by Reference 12 and it is formulated in Reference 1. Determining each 
resizing step under the ass\imption of constant modes (i.e., independent of the 
design variables), but using updated modes at each resizing step, may under- 
estimate the amount of material to be added locally for a certain amount of 
stiffening in a particular step, but it is not expected to cause an erratic or 
nonconvergent resizing path. 

As important as the frequency of updating is on the efficiency of the 
optimization procedure, the number of modes used to do the final flutter 
analysis is more important from a general point of view since it provides the 
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final check on the optimization procedure. In the opinion of the present 
investigators, a check flutter analysis using a proven sufficient number of 
vibration modes of the final configuration should conclude any optimization 
process. If flutter requirements are not met, then a new optimization process 
can be initiated and, probably, more modes or more frequent updating, or 
both, should be used. 


. 5 Ee commendations 


In view of experience diiring the present study, and as a result of 
experience with flutter analyses of actual airplane designs, the present 
investigators recommend the following: 

1. A flutter module should provide the option of inputting arbitrary initial 
modalizing matrices or of generating initial modalizing matrices based 

on a vibration analysis of the initial configuration. 

2. The number of vibration modes to be used for the flutter calculations 
should be an input option. 

3. The frequency of updating the vibration modes should be an input option. 

k. An option should be included to provide the analyst with information to 

determine whether his choice of number of modes and frequency of updating 
has led to satisfactory flutter characteristics. Such information might 
be provided by a vibration and flutter analysis of the final configura- 
tion with more modes than were used throughout the resizing process, 
a check on whether the optimality criteria for flutter are satisfied, 
or other check procedures. 


5 . AERODYNAMICS 
5 . 1 Introduction 


One of the objectives of this study is to develop general, efficient and 
accurate computational procedures for evaluating the vmsteady aerodynamic 
parameters necessary for use in a flutter optimization module, without, how- 
ever, evaluating aerodynamics theories. 

The procedure should be general. That is, it should be applicable to 
all present, and hopefully future, theoretical formulations of unsteady 
aerodynamics . 
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The procedure should be efficient. In the context of application in a 
flutter optimization module, this implies a minimum of computational operations 
required to recompute the generalized aerodynamic force coefficients each 
time a modal updating occurs. 

The procedure should be accurate. This implies it should be able to 
accommodate the most sophisticated formulations of the aerodynamics, such 
that the aerodynamics used in the flutter optimization module have the same 
accuracy as the aerodynamics used in a flutter analysis module. 

In the following section general background for a matrix formulation 
that allows a procedure satisfying these requirements is presented. It is 
followed by the definition of the formulation and a discussion of how the 
dimensions of the matrices, the method of interpolation for modal deflections 
and arbitrary values of the reduced frequency k, and the number of reduced 
frequency intervals to be considered determine the sequence of operations 
that is most efficient. Conclusions and recommendations regarding the aero- 
dynamics subroutine in a flutter optimization module are presented. 


5 . 2 General 


The elements of the matrix of generalized aerodynamic force coefficients 
are defined by: 


= J'J'f^(x,y) Pj(x,y) dx dy 


(5.1) 


Here p.(x,y) is the lifting pressure distribution associated with an 
angle-of-attack distribution, Q!,(x,y), which is defined by: 


«j(x,y) = 


ik df (x,y) 


9x 




(5.2) 


Zp 9z 

which expresses as the sum of — and terms in the case of harmonic 


dx 


motion with reduced frequency k in a mode defined by f.(x,y). 

J 

Expressed in the form of equation (5.l)»the evaluation of A. requires 

J 

evaluation of the surface integral each time new modes f^ are used. In the 

usual flutter investigation many different sets of modes are used, corre- 
sponding to different weight and stiffness distributions. In addition it is 
expected that frequent remodalization is required in an optimization procedure. 
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Thus it is advantageous to develop a method in which the generalized aerody- 
namic force coefficients are formed from a mode-independent part that contains 
as many of the numerical operations as possible, and a relatively simple 
mode-dependent part. 

Foiir different approaches are recognized in separating mode-independent 
operations from mode-dependent operations . One method relies entirely on 
analytical evaluation of the surface integral (Equation (5-1)). A second 
method formulates a numerical evaluation of the surface integral leading, 
effectively, to "lumped" aerodynamic forces at a grid of integration points. 

In a third method pressure distribution modes are analytically integrated 
over small areas and combined into elementary aerodynamic forces directly 
comparable to, and treated as, inertial forces. The fourth method recognized 
is based on a finite element approach, the basic formulation of which has no 
reference to pressiire distributions over the entire surface. 


The first three methods are usually thought of as stemming from the 
kernel function approach of Reference 13. In it the pressure distribution 


Pj(x,y) 
p’^(x,y) : 


is assumed to be a linear combination of pressure distribution modes 


P.(x,y) = Sa” p“(x,y) 
U J 


(5.3) 


The pressure mode coefficients 


n 


are determined from a boundary 


condition requiring that the normalized induced velocity distribution result- 
ing from the pressure distribution equals the angle-of-attack distribution at 
a set of downwash collocation points; 


K) ' 

Combining equations (5-3) and (5-^) leads to 

(pj) = (5-5) 

where the elements of matrix [PKI] are the integrals of the product of 
pressure distribution mode and an aerodynamic kernel. 

The columns of linearly independent pressure distribution 

modes . 
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5.2.1 Analytical Integration - When p.(x,y) is a linear comhination of 

J 

pressure distribution modes p^(x,y), analytical modal functions fj^(x,y) 
can he selected such that the integrals P"( x,y) dx dy can be 


evaluated analytically. Generalized aerodynamic force coefficients in terms 
of modal coordinates can then be formed. The analytical modes can be used 
as arbitrary modes to modallze the flutter equation, or a mode-dependent 
transformation between the analytical modes and the actual modes is used to 
express the generalized aerodynamic force, coefficients in terms of the actual 
modal coordinates . 


5.2.2 Numerical Integration of the Product of Displacement and Pressure — 
Reference 12 defines an approach to separating mode-independent operations 
from mode-dependent operations in which the surface integral of equation (5.1) 
is evaluated numerically. A Gaussian integration procedure is suggested to 


evaluate the integral. The pressure p.(x,y) and the deflection h.(x,y) 

are evaluated at integration points defined by the Gaussian procedure. 
Weighting factors in the form of a row matrix, [wfJ, make it possible to 
write : 


'^y ^ |^wfJ |h^ p^j 

The right hand side of equation (5.6) can be written as: 

M K = H h] ("J “ LhJ W hi 


(5.6) 


(5.7) 


The interchange of the row and diagonal matrix in the latter part of the 
combined equation (5.T) makes it possible to separate the mode-dependent 
operations from the mode-independent operations. 


The column matrix 


{wF.p.j 1 


is a set of lumped aerodynamic forces. The 


deflections h^ can be expressed in terms of the deflections z^ at the 
structural nodes by the relation. 


|h^(x,y)j = I^h] {z.j 


(5.8) 


A variation of this method is obtained if instead of the pressure 

distribution the velocity potential distribution <P. (x,y) due to an 

J 
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angle-of-attack distribution a^(x,y) is used. With the familiar linearized 
relation between pressure and velocity potential: 

5 = H < 5 . 9 ) 

equation (5.l) becomes: 


'll “ ^ 


dy 


(5.10) 


It can be shown that with the help of numerical techniques equation (5.10) 
can be written as: 


A. . 

ij 


hJ ["] ^ M] {'’j) 


(5.11) 


where [wf] performs numerically the first integration in equation ( 5 .IO) and 
[WFD]. performs the differentiation and integration in the second term of that 
equation. Equation (5.1l) is a triple matrix product, similar to equa- 
tion ( 5 . 7 ) j in which the center matrix is mode independent. For additional 
details see Reference 1. 

5 . 2.3 Numerical Integration of the Pressures - When p (x,y) is a linear 
~ J 

combination of pressure distribution modes p'^(x,y), the integral 

J^p^(x,y) dx dy can be evaluated over small areas, often referred to as 

aerodynamic boxes, into which the siirface is divided. By evaluating 

JJ ^ P^(x,y) dx dy and Sfy P*^(x,y) dx dy over the same areas, lumped 

aerodynamic forces can be determined in magnitude and location. The modal 
displacement at the location of each lumped force (i.e., for each aerodynamic 

box and each p'^(x,y)) can be expressed in terms of the structural degrees 

of freedom. Thus the product of each limiped force and its modal displacement 
can be formed. Summation over the aerodynamic boxes and the pressure dis- 
tribution modes participating in p.(x,y) leads to A . 

J i J 

5 . 2 .U Finite Element Approach - In a finite element approach, lumped aerody- 
namic forces corresponding to |wF*p^J (see equation (5.7)) are expressed 


2h 



under appropriate 


directly in terms of angle-of-attack distributions 



simplifying assumptions. The generalized aerodynamic force coefficients are 
formed as in Section 5.2.2. 


5.3 Basic Formulation 


Whatever the approach, or whatever aerodynamic theory is chosen, the 
generalized aerodynamic force coefficients in terms of modal coordinates 
can be expressed as the product of five matrices of which only the first and 
last are mode-dependent: 

hj] = H" [^r [""] ["] H 

The matrix [AIC] = [AlC(k)] , a function of the reduced frequency 

Ci)C 

k = -^ and the Mach number, is the core of the aerodynamics and is independent 
of mode shape. Its elements are basic aerodynamic influence coefficients 
defining lumped aerodynamic forces |^a.| aerodynamic force grid in 

terms of the angles of attack at downwash collocation points : 


K) = [AIC(k)]{cj 


(5.13) 


[w] = [w(k)J = [[dx] + ik [dz]] relates the angles of attack to the structural 

displacements {z} • It is independent of mode shape. 

IT 

is independent of k and of mode shape, and distributes 


The matrix 


[H]^ 


lumped aerodynamic forces and moments over the structural coordinates . 


In the case that the approach of Section 5.2.1 is followed, [AIC] is 
the matrix of generalized aerodynamic force coefficients in terms of the 
analytical modes; [H] and fwj are equated to 


r 1 

T r T~ 

^ r IT 


H = 

[f^(x,y)J [f£(x,y)J 

1 

(5.1U) 
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are 


where columns of 


[f£(x.y)] 


the fixed analytical modes. The matrix [f] 


is independent of k and of the actual mode shapes used to reduce the order 
of the flutter equation. 

The operations performed by [W] and may be included in [aic]. 

Equation (5.12) then reduces to the product of three matrices. 


The matrix [zj contains the modal columns in terms of the structiiral 
deflections {z}. 

The matrices [AIC], 

procedure. They will be used many times during the design process of an air- 
plane with a given external configiaration. It is therefore advantageous to 
form these matrices in a special aerodynamics computer program. 


and [w] are constant during an optimization 


Each time during an optimization procedure that a remodalization takes 

place, must be recomputed. Depending on the dimensions of the matrices 

in equation (5.12), it may be more efficient to compute the triple matrix 
product 


[haw] . [„]* [aic] [w] 


( 5 . 15 ) 

in the aerodynamics program, or to perform one or both of the multiplications 
HI? and [w][ z] in the optimization program. 

In the following sections factors are discussed that must be considered 


in determining which approach to nimierically evaluating 
equation (5.12) is most efficient. 


N 


according to 


5.^ Factors Affecting The Efficiency of The Numerical Evaluation of The 
Matrix of Generalized Aerodynamic Force Coefficients 


It is believed that the formulation 

H.] = [f H H 


( 5 . 16 ) 
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which follows from combining equations (5-12) and (5*15) » is widely used in 
industry. Detailed study of the formulation of equation (5»12), which 
directly follows from Reference 12, however, indicates that there are condi- 


tions under which it is more efficient to compute 


and/or [wjfz] 


in the flutter optimization module. Extensive comparisons have been made 
and are discussed in detail in Reference 1. In the following the factors 
affecting these comparisons are discussed and major conclusions are presented. 


5.^.1 Matrix Population - When matrices are sparsely populated, or populated 
in well defined blocks, proper programming can take advantage of this. 


In equation (5.12) the matrices [W] and .[h] may be sparsely populated. 
These matrices perform an interpolation and [w] , in addition, determines 
streamwise slopes at collocation points. In the Case of simple interpolation 
(linear or low degree polynomial), each row in [Wj expresses the angle of 
attack at a downwash collocation point in terms of several surrounding struc- 
tural coordinates . Similarly, each row of [H] expresses the deflection at 
an integration point in terms of several surrounding structural coordinates . 
Thus each row in [wJ and fnj contains relatively few, say <20, nonzero 
elements ; for linear interpolation each row contains four nonzero elements . 

In the case of interpolation by the surface spline method, the matrices [w] 
and [h] are fully populated, at least in the blocks that cover the aero- 
dynamic surfaces . 


Without specific stipulations, equation (5.16) implies that the order of 
the matrix [haw] equals the number of structural coordinates. There may 
be a considerable number of structural coordinates that do not carry an aero- 
dynamic load. They correspond to zero elements in [haw] . It is not expected 
that the fraction of nonzero elements will be high enough to justify treating 
[haw] as a sparse matrix. However, by proper ordering of the structural 
coordinates, the nonzero elements in [hAW] may be concentrated in one or 
more blocks . Then the aerodynamics program may form [hAW] based on aerody- 
namic load carrying structural coordinates only. Correspondingly, the flutter 
optimization module must eliminate structural coordinates that carry no aero- 
dynamic load from the modalizing matrix [z] . 

5.^.2 Interpolation for Arbitrary k Value - It is generally accepted that 
when the generalized aerodynamic force coefficients are determined for a 

discrete set of values, k^ , of the reduced frequency, inter- 

polation is adequate for approximating [ Aij(k) ] at arbitrary values of k. 

Two methods of interpolation are considered: cubic polynomial and cubic 

spline . 
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L 



can lead to "hunting" (oscillation between k values in adjacent intervals). 

by: 


Therefore it is recommended to define 

4 


where 




[a'jOo] 

[Ay(k>] 


( 5 . 20 ) 


t A!.(kfl) is the derivative of 

ij ^ J 

is an input to the interpolation subroutine. Thus the difference 


evaluated at . 


between equations ( 5 . 19 ) and (5.20) is that in equation (5.19) differentiation 
occurs after the polynomial fit and in equation (5.20) it occurs before the 


polynomial fit. The formation of 


^A^j(k^)J in the flutter optimization 

module is based on equation (5.12), equation ( 5 . 16 ) or any variant that is 

chosen as being most appropriate. The derivative [^AIC'(k)J or [HAW’(k)3 is 

needed and should be calculated outside the flutter optimization module by 
any method that gives adequate accuracy. 

To define A. .(k) and Al.(k)| in one k interval, foior matrices 

L J L J 

|^AIC(k^)J and four matrices |mC’ (k^)J , or four of each of the matrices 
[HAW(k)] and [HAW’(k)] must be input into the flutter optimization module. 

If [AlC(k)] and [AlC'(k)] are input, (k)J follows from equation (5.12). 

[a.j(k)] is given by: 

[A^j(k)j = i[z]^[H]^[AIC(k)J .(5.21) 

If k moves to an adjacent interval only two of the input matrices, one 
for the aerodynamics coefficients and one for their derivatives , need be 
replaced. 

It should be noted that, in effect, each k interval has its own 


associated polynomials for the value of 




and its derivative. 


The cubic spline method also defines different cubic polynomials for 
each k interval. The coefficients for the polynomial, however, are derived 
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£= li 2 . . . n. They follow 


from matrices defined for all k values k^ , 
from the assumption of continuity of derivatives over the complete range of 
k values. The resiiltant expression, e.g., for [AIC(k)J , is: 


[AIofj(k)]=[AICjo] + [AIOjJ(k-kp H- [AlCjjJk-k^^ ■. [AIC^3](k-kj 


r (5.22) 


The matrices 
optimization module. 
Then: 


[aic^„] to [aicJ 


should be formed outside the flutter 


[a|j (k)] = (k-k^) t (k-kp2 t |-AX^3j (k-kj)3 + 

|aZ^qJ ik + jAZ^^J ik(k-k^) + j^AZ^g] ik(k-k^)^ + I^AZ^J ^k(k-kg) 


(5.23) 


where 


and 


[“to] ' H’’H''["'^to]HH 

ho] = Hlf Ko]HH 


etc . , 


(5.2U) 


etc . , 


(5.25) 


Because of the implied continuity of the derivatives, it is proper to 
differentiate equation (5.23) directly and thus no additional matrices for 
the derivative need to be formed. 


To define 


namics input is 


hh 

is [aIC^(,] 


and 

to 


h‘"'] 


in one k interval if the aerody- 
requires eight coefficient matrices. 


Switching k to any other interval requires replacing all eight matrices. 


h'to] [“"m]’ 


If the basic aerodynamics input is in the form of 
then only four coefficient matrices are needed for each k interval. 

5 .U .3 Number of k Intervals - Let the basic aerodynamics input into the 


flutter optimization module be 
considered is 2 . 


[HAW(k^)j; 


the nimiber of k intervals to be 
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For the cubic polynomial, Lagrange's interpolation formula is considered 


to be most efficient since it expresses 


value 






directly in terms of its 


at discrete values k/- 


£ = 1 , 2 , 3 , k: 


(5.17) 


where <£^(k) is defined by: 


p ^ ^ (k-k^) (k-k^) (k-kj^) 

^~T(k) ~ 737 V 3 737 "37 S 737 _v 


(ki-k^) (k^-k^) (k^-kj^) 


(5.18) 


and cyclic substitution leads to '£)|‘ 

The interpolation formula (5.1?) is used only for the interval kg< k<k^. 

For the interval k^< k < kj^ the index 2 jnust be increased by one and for 

k^<k<kg the index must be lowered by one. 

Since most methods of optimization require the computation of the deriv- 
ative of the aerodynamics matrix with respect to k, the formation of the 


derivative 


’ [AijO')]. 


must be considered. 


Differentiating equation (5.17) with respect to k leads to an 
expression: 




(5.19) 


that is based on the same aerodynamic matrices as equation (5.17). This 
approach, however, combined with the re-indexing of k^ as k moves to an 
adjacent interval, leads to jumps in the value of K .(k)] at all values 
k=k^ . Apart from considerations of accuracy, this is undesirable since it 
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and 


If cutic polynomial interpolation is used, (2+3) matrices [HAW(k^)J 
|^HAW'(k^)J 


must be input and pre- and postmultiplied by 


[f 


and 
needed in 


l^zj to form the 2(£+3) matrices 
Lagrange's interpolation formula. 

If cubic spline interpolation is used, 4£ matrices to 

must be input and pre- and postmultiplied by form 


coefficient matrices 


ho] [^'43} 


Under othen/ise equal circumstances, polynomial interpolation is more 
efficient if 


2(£+3) < hi or 2>3 ( 5 - 26 ) 

This condition is valid for other sequences of operations to form 

(k)J according to equation (5-12) . However, there are also sequences 

for which > 4 or 2 > 5 is required for cubic polynomial interpolation 
to be more efficient than cubic spline interpolation. 

The number of intervals that should be used is difficult to predict. 

If only one flutter constraint is active, k may stay within a rather small 
range during the entire optimization process and that range may lie completely 
within one k interval. Obviously only aerodynamic matrices applicable to 
that one k interval need be computed. In general, however, several k 
intervals are required. 

5 . 4.4 Sequence of Multiplications - Defining one computational operation as 
one multiplication and one addition , the numbers of such operations required 
inside the flutter optimization module for different sequences of multiplica- 
tions in equation (5.12) have been determined and compared. 

The following options have been considered; the numerals Indicate the 
sequence of multiplication. 
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HI: 


H2; 


M = •["]"■ H-H'H 

' V ' ' V ' 

1 2 

. ' 

.3 

V / 

h 

« ' > / 

1 2 

' V ' 

3 


1* 


(5.27) 


(5.28) 


H3: 


[AIC] • [w] = [aw] is computed outside the flutter optimization 
module 



(5.29) 


Hl+: 



is computed outside the flutter module 


1 2 

> ' 

3 


(5.30) 


H5: 



is computed outside the flutter module 


KJ - H" H H 


(5.31) 


The number of computational operations is independent of sequence of 
multiplications in H5. 
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Formulas defining the number of numerical operations have been derived 
and are reported in Reference 1. Wo option stands out as clearly superior 
or inferior to all others, but some comments are offered in Section 5*5. 

5 . . 5 Fom of Inputting the Angle-of-Attack Generating Matrix - In the pre- 
vious section the options are defined as if one complex matrix [W(k)] was 
input for each value of k. Suboptions of the options that require inputting 
[W] can be obtained by considering the definition of [w]. 


[W(k)] = [[dx] + ik [dz]] 


(5.32) 


Thus in options HI-, H2 and Hi; an "a" and a "b" version can be recognized. 
In the "a" version [w(k)] is input for several k values. In the "b" 
version the real matrices [dxJ and [dzJ are input. 

Whether option "a" is more efficient than option "b" depends on the 
other factors. In the case of option HI it seems that Hla is favored if 
interpolation for deflections is held simple. Option Hlb is favored if 
surface spline interpolation is used. If the "a" option is used and the 
derivative of the aerodynamics matrix is needed, matrix [dz] must be 
input anyway. 


5.5 Summary of Comparisons 

Detailed comparison of the options HI through H5 is reported in Refer- 
ence 1. The following summarizes the comparisons. 

5 . 5.1 Input Storage Requirements - If the number of k intervals to be used 
is three or more, cubic polynomial interpolation for arbitrary values of k 
requires less input storage than cubic spline interpolation for all options 
HI through H5. 

5 . 5.2 Core Space - For options H5 and H3, cubic polynomial interpolation 
for k requires two times as much core space as cubic spline interpolation. 
For all other options, both methods of interpolation require the same core 
space . 

5 . 5.3 Read-In - Cubic polynomial interpolation for arbitrary k requires 
less read-in than cubic spline interpolation as the value of k moves into 
an adjacent interval. 
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5.5.^ Uumber of Computational Operations - Cubic polynomial interpolation 
for arbitrary values of k reqiiires fewer computational operations than cubic 
spline interpolation under the following conditions: 

options H3 and H5: if the nxmber of k intervals is more than three. 

options HI and HUi if the nxjmber of k intervals is more than fom:. 

option H2: if the number of k intervals is more than five. 

VJhich of the options HI through H5 is most efficient depends strongly on 
the dimensions of the matrices. These, in txirn, depend on the desired accu- 
racy, and the methods used for integrating or lumping the aerodynamic presstires 
and interpolating and differentiating the structural displacements . 


Let M be the number of modes, W the number of structiiral coordinates 
before modalizing, D the mrniber of downwash collocation points and K the 
nxmiber of integration points. Then equation (5*12) can be annotated as follows 


M = [f [*7 M H M 

M,W N,K K,D D,N N,M 


(5.33) 


with 

and 
to N 



From this equation it can be seen that if K and D are small compared 

H it becomes advantageous to perform the multiplications ’ 

[wj • [z] in the flutter optimization module. If K and D are equal 
or larger, then it becomes advantageous to form the product 

[aic] [w] outside the flutter module. The relationships defining 



when an option is better than another are complicated. They are documented 
in detail in Reference 1, but they have not led to simple criteria. 


5.6 Conclusions and Recommendations 


The preceding sections lead to the following conclusions: 

1. Formulation of the generalized aerodynamic forces in the form 

[Ay(k>] = [s]’'[H]’'[AIC(k)][w(k)][i] 

is possible and practicable for all approaches to determining unsteady 
aerodynamic forces. 


2 . Several options of inputting the matrices [h]'^ , [AlC(k)] , and 

[W(k)] = [dx] + ik[DZ] can he recognized, i.e., separately, after multi- 
plication, in the form of cubic spline matrices, or in derivative form. 
Which option is most efficient depends strongly on sizes of the matrices, 
whether the derivative of the aerodynamics matrix is needed, the method 
of Interpolation for arbitrary k, the population of the matrices [H] 
and fw] and the number of k intervals expected to be active during 
the optimization process. 


3. 


In addition to depending on the number of computational operations 


required to form 


[a.j(R)] 


the efficiency of the flutter optimization 


process may also depend on the required input storage, the possibility 
of storing all matrices required for interpolation of the aerodynamics 
for k in one interval in core, and the read-in required if the value 
of k moves to an adjacent interval. 


4. It is possible to design a flutter optimization module that includes all 
options such that the user can choose the option that is most efficient, 
that fits his available data or that he prefers for some other reason. 


In view of the above conclusions the following recommendation is made: 


In designing a flutter optimization module for a facility, the calculation 
of the generalized aerodynamic force coefficients and their derivatives should 
be based on the formulations presented in this section. That is, a mode- 
independent part should be generated outside the flutter optimization module 
leaving an often to be repeated mode-dependent part of the calculations to be 
performed inside the flutter module. The number of options is large and it 
may not be practicable to include all options in the flutter module. The 
choice of options is facility dependent. Certain practices of generating 
generalized aerodynamic force coefficients may already exist and the existing 
computer system may influence the choice significantly. The module, however, 
should allow the user considerable freedom in choosing the option that is most 
efficient for his problem. The number of options to be included should be 
decided on the basis of a stand-alone flutter optimization module. It should 
not be restricted because of the module being part of a general analysis system 
which at present has only a restricted choice of outputting aerodynamic 
coefficients. 
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6. METHODS OF OPTIMIZATION FOR FLUTTER 


6 . 1 General 

In this section five approaches to structural optimization with flutter 
constraint are reviewed. All five belong to the category- of direct methods'; 
i.e. methods in which the mathematical formulation defines resizing steps 
aimed directly at determining the extreme value of the objective function, 
in this case minimizing the structural mass. In contrast, in the indirect 
methods the mathematical formulation defines resizing steps aimed at satisfy- 
ing a criterion that, when satisfied, indicates that the optimum condition 
is reached. 

The five methods reviewed represent distinctly different mathematical 
formulations. They are, in the order of review: the gradient methods of 

Rudisill and Bhatia (Reference l4), the weight gradient method of Simodynes 
(Reference 15), a penalty function method (Reference l6), a method of 
feasible directions (Reference IT), and a method that evolved from the method 
of Incremented Flutter Analysis (References and l8). All are formulated 
under the assumption of modalization matrices that are independent of the 
design variables, but can be updated at any resizing step. 

In addition to the qualitative evaluation of the methods presented in 
this section there is a numerical evaluation and comparison in Appendix A. 

In that appendix the results are presented of applying these methods to an 
optimization task in which the flutter equation is written in terms of h9 
discrete degrees of freedom and is not modalized. 


6.2 Rudisill-Bhatia Approach 

Reference l4 defines four different resizing columns that alone or in 
combination can be used to design a minimvim weight, or near minimum weight, 
structure that satisfies a flutter speed constraint. The resizing columns 

are defined in terms of incremental values, "the design variables 

. (The notation of Reference l^t is followed.) Three of them involve the 

gradient of the flutter speed with respect to the design variables. In con- 
nection with this. Reference lU presents closed form analytical expressions 
for the derivative of the flutter speed 'vrilth respect to a design variable. 
These expressions have been used successfully in n'umerical test cases per- 
formed during this study. N'umerical values of the derivatives are in good 
agreement with values obtained using the approach of Reference 15. In the 
following the resizing columns defined in Reference iH are discussed. The 
terminology of Reference li+ is followed. 



6.2.1 Velocity Gradient Search - Equation (23) of Reference lU defines a 
coliomn of design variable increments as: 


AP. 

1 

AV 

av 

n /_ \2 

ap. 



1 1 


( 6 . 1 ) 


8V 

Equation (6.1) defines a direction by means of 


The magnitude 


of the increments is directly related to the velocity increment AV through 

AV 


the coefficient 


k fty 


, where AV = 


= ^ I {ap 

KJi ' 


is the approximate 


increase in flutter speed due to the design variable increments 


{^"i) 


The formulation is known as the method of steepest ascent of the function 


V(P.). 

cosines 


The direction implied by equation (6.l) 

f 'I f 'i 


dP. ■ 
1 

V' "uVl "1 V\ 1 r ” 

dV 

dS 

tor which - 

) 


defines a column of directional 
is maximum. 


Assuming a linear relationship between the design variables P^ and the 


associated mass , 


m. 

1 


"i"i’ 


equation (6.l) can be written as: 



( 6 . 2 ) 


Taking ^ ^ as a reference it can be seen that the direction of 

depends on the scaling between the design variable and the associated mass . 
Since the quantity of direct interest is total mass, and not design variables 
related to mass, it seems logical to choose elementary mass as design 
variables and choose that case equation (6.1) becomes: 



which represents a design change in the direction of maximtmi 


(6.3) 
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The usefulness of the velocity gradient search is related to its ability 
to raise the flutter speed more efficiently, i.e., with a smaller increase in 

total mass m^ than a simple increase of the overall stiffness level. Simple 

physical considerations, however, lead to the conclusion that the most 

dV dV 

efficient move is in the direction in which ^ is maximum. 

dM rdm. 

1 

The nimierical evaluations in Appendix A indicate that a resizing colimin 


proportional to |9m. | efficient means of resizing a structure in one 

step to satisfy a flutter speed constraint with a moderate mass penalty, thus 
providing a good starting point for a procediure that minimizes the total mass 
at constant flutter speed. 


The relative efficiency of a resizing column proportional to |9m.| 

follows from the fact that it tends to add more material where it is most 

8V 

efficient in raising the flutter speed. Design variables for which -5— is 

dm . 

1 

negative are reduced in value, which also raises the flutter speed, 

6.2.2 Mass Gradient Search - Equation (30) of Reference l^i defines a column 
of design variable increments : 


i ^ |9M 

J 9P. 8P. 

L iJ I 1 


(,6.U) 


where M is the total mass: M=Zm. . 

1 r ^ 

Equation (6.4) defines a direction by means of . The magnitude of 

i 

the increments is directly related to the velocity increment AV through the 


coefficient 




8p. lap. 
L iJ t iJ 


The direction defined by equation (6.4) corresponds to a maximum value of 


'M' 
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Again assuming a linear relation between the design variables and 


the associated mass m.=C.P. it can be seen that 

1 XI 


n 


n 


M = 2 = 2 C.P. 

1=1 1=1 


(6.5) 


and ^ - C. . 


Thus equation (6.U) becomes; 


K) = 


AV 

c- F- K1 

1 9m. 1 

L ij ^ J 


Kl = 


AV 


V p2 8V 

."^i 8m. 
1=1 1 


K) 


( 6 . 6 ) 


It is apparent that the direction of again depends on the scaling 

between design variable and associated mass. 

Choosing m. as design variables equation (6.4) becomes; 


W 


^ 8V 

.A, 8m. 

1=1 1 


(6.7) 


which represents a uniformly distributed weight increment. 

Reference l4 us'es equation (6.4) to define a decrease in total mass, 
thus a negative AV is used. By following a path of steepest descent for 

M(P^) the emphasis is on decreasing total mass. The relation AM vs. AV 

is not considered. 

The mass gradient search could be used in combination with the velocity 
gradient search to formulate an optimization procedure. Alternate application 
of these searches tends to lower the design variable weight required for 
satisfying the flutter constraint since the velocity gradient search tends to 
Increase the flutter speed by adding a relatively efficient mass distribution. 
Although the mass gradient search removes mass indiscriminately, repeated 
application of both searches tends to an optimum mass distribution. Rather 
than remove mass indiscriminately, it seems logical to remove mass first where 
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it reduces flutter speed the least. That would be the case if most mass would 
9V 

be removed where is smallest. Note that if there are design variables 

9V ^ 

for which is negative their reduction increases the flutter speed, 

i 

allowing more mass to be removed from design variables with a small positive 
9V 

value of . An efficient method to remove mass is to remove it propor- 

^1 9V . 

tionally to 9 v/^m '' design variables for which is positive. 

In Reference 1 the sequential application of a velocity gradient search 
and a mass gradient search is replaced by the application of an equivalent 
two-component column of design variable increments. Regrouping of the ele- 
ments in the two components shows clearly that this method does result in a 

9V 

resizing where more mass is added where - — is larger and more mass is 

dm. 

9V ^ 8V 

removed where ^ — is smaller, regardless of the algebraic sign of r — . 

dm. dm. 

1 1 

6.2.3 Gradient Projection Searches - The gradient projection search follows 
a direction of steepest ascent while satisfying a constraint. The gradient 
projection search used in Reference l4 is aimed at following the steepest 
ascent of the flutter speed as a function of design variables while keeping 
the total weight constant . 

From equations (32) and (3^+) in Reference l4 the following column of 
design variable increments can be derived: 


AP. 

1 


AS 


\/y \^+ -.mm 

1 9P. 9P. 

\ 1 / 11 


9P. 1 9P. 

1 1 


( 6 . 8 ) 


2 2 

where AS = 2 ('^P- ) > the step size in terms of design variables, 

i=l ^ 

Equation (6.8) can also be written as 


AP. 

1 

AV 



9V 

9P. 

. i_ 


9V,., 9M 1 

9P. ^1 9P. 

1 iJ 

9P. ‘'■1 9P. 


Uo 


( 6 . 9 ) 



In "both equations 



which follows from the condition 



( 6 . 10 ) 


( 6 . 11 ) 


Usually V(P^) has a maximum for a given total weight. Although in the 

application shown in Reference l4 this maximum is not sought, it may be 
reached. It would seem, therefore, that the formulation of equation (6.8) with 
AS detemiining the step size, is preferable over the formulation of equa- 
tion ( 6 . 9 ) in which AV determines the step size. 

Reference l4 mentions the possibility of a gradient projection search 
in which the flutter speed is held constant and the total weight is reduced. 

The present authors consider this a more significant procedure from a prac- 
tical point of view. As indicated in Reference l4 the formulas for this 
approach can be obtained from the constant total weight approach by inter- 
changing the symbols V and M and changing the algebraic sign on AS. 


Since the total mass M has a minimum value, only the equivalent of 
equation (6.8) will be given: 


where : 



9M av 

8P. ^1 9P. 

1 1 


9M ^ 

9P. 1 

■V _ 4. ij 

9V ' 

ap. 

i ij 

'^'l 

9V 

9P. 

L iJ 

[avl 

18P.^ 


which follows from the condition 



(6.12) 


(6.13) 


(6.14) 

4l 



Assioming a linear relationship between the elementary masses, m^, and 
the design variables , : 


m. = C.P. 
1 11 


(6.15) 


equations (6.12) and (6.13) become: 

-AS 


AP. 




1 


(6.16) 


X = 

1 


L 1 J 


dP. 

1 


dp. 

1 


dV ' 

dP. 

1 


(6.17) 


Letting, in addition, = 1, i.e., choosing the elementary masses as 


design variables, these two equations become: 


Am. 

-AS 

i*x, 

1 

N /n + X, S ' 

1 dm^ 


1=1 1 


( 6 . 18 ) 


^1 = 


n 

2 ^ 

i=ldmi 


n 

i?iV 


/il' 

dm. 


(6.19) 


The step size is determined by 


P ^ P 

IAS = 2 Am. 

' itl ^ 


( 6 . 20 ) 


h2 


r 


Equation ( 6 . 18 ) will now be discussed under the assumption that 

2 > 0 and thus < 0 . This implies that j^^.J |^| latter 

inequality indicates that a uniform increment of all the design variables 
increases the flutter speed. This is not an tinreasonable ass\miption. Addi- 

” 8V 

tionally it can be shown that n + X^ " ^ ~ arbitrary set of 

8V ^ 

values of . 

i 


In view of these considerations equation ,( 6 . 18 .) can be written as: 


where = 


than 1. 



r ’] 






Am. 

1 

= 

-^1 


1 

AS 


>0 

and 


8V 

1 

f n 

. *v ^ 

8V 

^1 


( 6 . 21 ) 


1 3 m. 

1=1 1 


The resizing coliunn defined by equation ( 6 . 21 ) is a linear combination of 
the columns associated with the velocity and mass gradient search discussed in 
the previous sections. The algebraic sum of the two contributions tends to 
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increase the design variables with a larger value of 
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those with a smaller value of , including those with a negative value of 

8V ^ 

5— . This agrees with physical reasoning. As defined in equation (6.I8) the 
two contributions to equation ( 6 . 21 ) are in a fixed ratio. Due to 

nonlinearities this results in a drift in flutter speed due to a resizing 
step • 
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an equivalent of equation (6.21) is formulated that allows the use of the 
method of_Incremented Flutter Analysis (References and l) for determining a 
value of K such that the velocity increment, AV, associated with the velocity 
gradient column is exactly cancelled. 
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In equation (6.21), determines the step size. In a follow-on paper 
(Reference 19) Rudisill and Bhatia describe a deterministic method of 
defining step size in the case that the gradient projection search is used 
to increase the flutter speed at constant total mass. The method is based on 
the assumption that the flutter speed is a nearly quadratic function of the 
design variables . 

Since the total mass is a linear function of the design variables this 
method of determining the step size is not applicable to the gradient pro- 
jection search at constant flutter speed. Choosing a step size then is a 
matter of experience and judgment. 

6.2.U Concluding Remarks - The preceding discussion deals with the basic 
resizing columns as defined by foirr different approaches. In any practical 
application, minimum size constraints must be taken into account. Conceptually 
this is a simple matter. Computer programming requires close attention to 
detail since a variety of potential minimum size constraints, some of them 
expected to occur infrequently , must be foreseen. Reference l4 does not go 
into detail on this . 

The capability of the approaches presented in Reference li+ for structural 
optimization with a flutter constraint lies exclusively in the iise of the 

of the flutter speed with respect to the design variables . Repeated 
of structural material according to a distribution proportional to 

followed each time by a rather indiscriminate removal of material to 

the desired flutter speed, converges to a minimum mass design. 
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6.3 The Weight Gradient Method of Slmodynes 


In March of 1973, E. E. Simodynes presented a method for the optimiza- 
tion of structural weight at a specified flutter speed (Reference 15). The 
method employs the gradient of total weight with respect to n-2 design 
variables in a resizing algorithm to minimize structural weight while main- 
taining a constant flutter speed. Two of the design variables are dependent 
variables. The resizing column is such that a constant flutter frequency is 
also maintained, a characteristic which simplifies the formation of the 
required derivatives. The restriction on flutter frequency, however, repre- 
sents an arbitrary constraint on the procedure which may result in weight 
penalties which can not be justified by the computational simplifications 
achieved. A modification of the method is proposed wherein the flutter 
frequency is permitted to vary, taking the place of one of the dependent 
design variables. In the following, the original method is discussed first, 
followed by a discussion of the modified procedure. 

6.3.1 The Method of Simodynes - The optimization method presented in Refer- 
ence 15 uses a resizing algorithm based on the gradient of total weight. 


1 + 1 + 
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subject to flutter speed and flutter frequency constraints . As with most 
optimization methods, the stiffness and inertia terms are expressed in 
linear form: 


The matrices 


[k] = [kJ + 2 and = [mJ + £ [.^mJ ( 6 . 23 ) 
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flutter equation for neutral stability is 


represent the stiffness and inertia of the 
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design variable. The 
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where [Q] is a matrix of unsteady aerodynamic force coefficients . 


Considering flutter frequency and flutter velocity to be specified, the 
fixed and variable coefficients of the flutter equation may be grouped as in 
equations (6.25) and (6.26), and the flutter equation expressed as in 
equation ( 6 . 27 ). 


[a] = [k„] - C 0 ^[mJ - 

(6.25) 
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Two of the design variables are now designated as the dependent design 

variables m and m , and derivatives of the flutter equation with respect 
u V 

to the independent design variables m^ formed as shown in equation ( 6 . 28 ). 
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In this equation, -[q} is the flutter eigenvector of equation (6.27), 
and j_rj is the corresponding row eigenvector. Since the complex coeffi- 
cients of equation ( 6 . 28 ) are known, the partial derivatives 9m^/9m^ and 
9m^/9m^ are readily obtained. The weight of the structure is now e3q>ressed 
in equation ( 6 . 29 ), where is the weight of the fixed structure; the 

derivative of the weight with respect to each independent design variable 
is computed (equation (6.30)), and the gradient of total weight formed 
( equation (6.31)) . 
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A resizing column of increments is then formed for a specific weight reduc- 
tion, -AW (equation (6.32)), and new values of m^ and m^ are calculated 
using this column of increments and the partial derivatives 9m^/9m^ and 
9m /9m. . 

V 1 



(6.32) 


6.3.2 Discussion of the Method - As indicated in Reference 15, the principal 
feature of the method of Simodynes is the weight gradient, and the resizing 
column derived from it. The computation of this weight gradient is partic- 
ularly straightforward due to the constraints on flutter speed and flutter 
frequency imposed on the resizing procedure. In particular, the unsteady 
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aerodyn^ic parameters are constant throughout the optimization cycle, so 
that the flutter derivatives do not involve derivatives of the aerodynamic 
parameters . 

This simplification, however, is not obtained without an associated 
disadvantage. The imposition of the frequency constraint on the optimization 
process results in an artificial constraint on the distribution of the design 
variable masses during the resizing procedure. An example of this may be 
seen in the numerical evaluations reported in Appendix A. The effect of 
this constraint is difficult to assess, since it depends to a great degree 
on the particular situation. For a configuration far from the final optimum 
distribution, with a frequency of the critical mode significantly different 
from the final value, the effect of the frequency constraint would be expected 
to be significant. For a configuration closer to the optimum distribution, 
with a critical mode frequency approximately equal to that of the final con- 
figuration, the effect would be correspondingly less. Reference 15 states 
that a fixed set of vibration modes is used throughout a complete optimization 
cycle, and it is concluded that the flutter frequency constraint is maintained 
throughout the optimization cycle also. In that case, at least one additional 
optimization cycle must be performed in order to assure that no effects of the 
constraint remain. The choice of dependent design variables also influences 
the magnitude of the constraint effect, since design variables having little 
influence on frequency would produce relatively greater distortions of the 
unconstrained distributions. 

As with other methods employing an arbitrary step-size (e.g. Refer- 
ence 1^), this parameter must be established on the basis of judgment, 
intuition, experience or (probably) a combination of these. For the present 
method, step-size is determined by the weight reduction, -AW, specified 
by the user. The choice of this parameter involves a compromise; too large 
a value of AW could result in unacceptably large exc\irsions in flutter speed 
and frequency; too small a value of AW would result in an excessive number 
of resizing steps required to reach optimum. 

Minimum size constraints are handled in a straightforward manner; when 
a design variable reaches minimum size, it is temporarily eliminated as a 
design variable. In subsequent steps, the derivatives are calculated for this 
design variable and it is reinstated as an active design variable if these 
calculations so indicate. Although not specifically stated in the reference, 
this is presumably true of the dependent design variables as well as the inde- 
pendent design variables. Here again, a poor choice of dependent design 
variables , requiring frequent shifting to other design variables during the 
optimization cycle, would result in an inefficient resizing process. 

6 . 3.3 A Modification of the Method of Simodynes - As discussed earlier, the 
constraint on the flutter frequency provides a simplification in the computa- 
tion of the required derivatives, but results in a nonoptimum weight incre- 
ment implying a weight penalty that is difficiilt to predict. A modification 
is suggested which eliminates this frequency constraint. This is done by 
allowing the flutter frequency to become a dependent variable, taking the 



place of one of the dependent design variables. Following a procedure 
similar to that of the original method, derivatives of the flutter equation 
are obtained after first expressing the flutter equation as a function of the 
reduced frequency, k. The result is equation (6.33), which is comparable to 
equation (6.28) of the original procedure. Equation (6.33) is solved for the 
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unknown parameters 9m^/9m^ and 9k/9m^; the total weight is expressed in 
equation (6.34), and the derivative of the weight with respect to each 
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of the independent design variables is shown in equation (6.35). 

The gradient of the weight and the colijinn of increments are formed as be- 
fore (equations (6.31) and (6.32)). Except for this modification in the weight 
gradient formation, the optimization cycle proceeds as in the original method. 

6.3.4 Discussion of the Modified Method - The modification of Simodynes' 
method suggested here unquestionably results in a less approximate procedure. 
The frequency constraint, with the associated distortion of the resulting 
mass distribution, is removed. Not only does this eliminate an iindesirable 
feature of the original method, but fewer optimization cycles should be 
required in order to reach a satisfactory approximation of the optimum dis- 
tribution. As a result of the modification, however, some additional complica- 
tion of the computational procedure is required. In solving the flutter 
equation to obtain the eigenvectors, a fixed matrix of unsteady aerodynamic 
coefficients can no longer be used throughout the optimization cycle. In 
addition, the derivatives of the aerodynamic parameters with respect to the 
reduced frequency, k, must be obtained. How troublesome these complications 



are depends on the versatility of the computational system used. For the 
numerical evaluation of the modified method presented in Appendix A, the 
Lockheed p-k method of solving the flutter equation was used; this procedure 
has a built-in subroutine for interpolating matrices of aerodynamic coeffi- 
cients . A similar routine was used to obtain the approximate derivatives of 
the aerodynamic parameters using a finite-difference procedure. Since these 
computational tools were readily available, the modifications resulted in no 
peirticular computational diffictilty. 


It should be noted that a certain similarity exists between the column 
of increments resulting from the present procedure and the, column of incre- 
ments resulting from the gradient projection search with constant flutter 
speed of Rudisill-Bhatia (Section 6.2.3). The expression for the column of 
increments of the present procedure is shown in equation (6.36) and from 
equation (A. 6) of Appendix A it is seen that equation (6.37) is an equivalent 
expression. 
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Recognizing that the scalar premultiplier of the column is arbitrarily 
chosen, the increments for the independent design variables can be made 

identical to those of equation (6.21) of Section 6.2.3 if the normalization 
factor l/aV/8m^ is equal to of equation (6.21) . It will be recalled 

that is chosen so as to result in a flutter velocity increment, on a 

linear basis, equal to zero. In general, the normalization factor for the 
present procedure will not be equal to since it results from the choice 

of the dependent design variable. The (linear) flutter speed increment is 
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then held to zero by the increment in the dependent design variable defined 
in equation ( 6 . 37 ). 


6.3.5 Assessment of tllb Method - The original Simodynes method is a simple, 
straightforward procedure which involves a minimum of computational difficulty 
in the generation of the flutter derivatives . The flutter frequency con- 
straint might be a 'useful device in adopting existing routines for the solu- 
tion of the flutter equation and generation of unsteady aerodynamic parameters 
for use in a flutter optimization procedure. In the development of an inte- 
grated design procedure, however, there would appear to be no clear advantage 
in retaining the frequency constraint. Instead, a modification of the pro- 
cedure such as indicated here would be highly desirable. 


S.h An Interior Penalty Function Method 


In seeking to evaluate a penalty function method, the initial task is 
one of defining both the method itself and the scope of the evaluation. As 
used here, the term penalty function method refers to any structural optimiza- 
tion technique in which penalty terms, which are related to the constraint 
equations, are added to an objective function to form a modified objective 
function, which is then minimized. For the purposes of this evaluation, 
however, only a single, representative procedure will be considered. The 
particular procedure chosen is an interior penalty function method and follows 
closely the method described in Reference 16, and is essentially identical to 
that used in the numerical evaluations presented in Appendix A. A brief 
description of this method is presented in the next section, preparatoiy to 
the discussion and evaluation which follow. 


G.h.l Description of Method - As indicated previously, the method considered 
here is based on the method described in Reference I6. In that procedure, a 

modified objective function, P(m^,r), is formed as shown in equation ( 6 ; 38 ). 


P(m^,r) 


W(m.) + r *2 — ^ 
1 £ g^(m.) 


( 6 . 38 ) 


The term W(m^) is the quantity to be minimized, or objective function, 

and represents the weight associated with the design variables m^ . The 

second term of equation (6.38) expresses the penalty terms as functions of the 
design constraints; r is a penalty function weighting factor. By means of 

this formulation, the problem of minimizing W(m^), subject to the constraints 

g^(m^), is transformed to one of an unconstrained minimization of P(m^,r) 
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using a SUMT (Sequence of Unconstrained Minimization Technique) approach 
(Reference 20 and 2l) . One such minimization is carried out for each succes- 
sive reduction of the value of the penalty function weighting factor, r, 

until the minimized value of P(m^,r) is approximately equal to the corre- 
sponding value of W(m^) . 

The unconstrained minimization of P(m^,r) is accomplished by first 

generating a move-vector direction, then determining the minimum of P(m^,r) 

in this direction by means of a one-dimensional search. For determining 
this move direction, several direction-generating algorithms are available, 
the best known of which is the DFP algorithm (Reference 22) as used in the 
preliminary design procedure reported in Reference 9- Reference l6, how- 
ever, uses a variation of Newton's method wherein the second derivatives of 


P(m^,r) are approximated, and this procedirre will be used here. The second 
derivatives shown in equation (6.39) may be approximated by neglecting 
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the second term of the summation, under the assumptions stated in Reference 

l6. In addition, the objective function, W(m^), is assumed to be a linear 

function of the design variables, so that the first term of equation (6.39) 
is equal to zero. The second derivatives are then approximated as in equa- 
tion (6.U0), and an estimate for the colimin of design variable 


aS 

am.am. 

1 J 


G. 

1 


= 2r 



(6.40) 


increments which will minimize 


P(m^ ,r) 


is formed as shown in equation (6.4l). 



(6.4i) 
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It should be noted that equation (6.4l) would exactly define the required 

column of design variable increments if P(m^,r) were a quadratic function 

of the design variables and the exact second derivatives used. In the present 
case, however, equation (6.i+l) is used only to define the move direction for 
the one-dimensional search. 

6.h .2 Discussion of the Method - In evaluating the penalty function method 
described here, four principal elements or characteristics can be discerned: 
the treatment of the constraints, the penalty term weighting factors, the 
direction generating algorithm and the one-dimensional search. In the fol- 
lowing sections, each of these characteristics will be discussed, followed 
by an assessment of the method as a complete procedure. 

The treatment of constraints is the principal distinguishing feat\ore 
of the penalty function method, and this treatment provides a number of 
advantages. As used in the present procedure, the inequality constraints 
serve two important functions: l) conditioning of a move vector such that it 

tends to avoid violation of the constraints, and 2) limiting of the move 
vector amplitude such that the constraints are not violated. The first of 
these functions is implemented through the direction generating procedure, 
and the second of these is a part of the one-dimensional search. 

The treatment of the constraints as inequality constraints, resulting 
in the characteristics described above, can be a very powerful approach in 
structural optimization for flutter. One of the more obvious advantages is 
the fact that the inclusion of multiple flutter speed constraints causes 
little conceptual difficulty. Whether this advantage is of practical value 
is difficult to assess without further work. To include all flutter speed 
constraints (for several Mach numbers and airplane loading conditions) in 
the penalty term would be a large computational burden. Thus it would seem 
that only active constraints should be included in the penalty term and that 
separate program logic should be used to determine which constraints are 
active. Experience at the Lockheed-California Company, partly obtained during 
this study, indicates that multiple active flutter speed constraints may occur 
rarely. Probably the most important advantage resulting from this method of 
handling constraints is the fact that a large number and variety of constraints 
can readily be included in an automated procedure. A disadvantage, however, 
is the fact that derivatives of the constraint quantities with respect to the 
design variables must be obtained. 

The handling of the penalty term weighting factors indicated in equa- 
tion ( 6 . 38 ) exerts a significant influence on the resultant performance of 
the method. The treatment of these factors is a matter of Judgment and 
depends on experience. Experience at NASA, Langley Research Center suggests 
an initial value for these weighting factors which makes the penalty terms 
approximately equal to the value of the objective function (Referenae 23). 

The final value for these factors can be determined by establishing acceptable 
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values of the residual constraint inequalities and then specifying an allowable 
percentage of the modified objective function contributed by the penalty terms. 
Having thus established the range of the weighting factors, the appropriate 
reduction factor to be applied during each step is determined from the selec- 
tion of the desired number of steps. Monitoring of the progress of the weight 
minimization, however, may lead the analyst to interrupt the optimization 
procedure and modify the penalty weighting factors. Specifically, optimiza- 
tion steps may continue after the number of steps on which the reduction 
factor is based are completed if the total weight progression indicates that 
a minimum weight has not been attained. 

At each step in the optimization process, the weighting factors deter- 
mine the extent to which the constraints influence the resultant move. If 
the weighting factors are large relative to the distance from the constraint, 
the design tends to move away from the constraint, into the feasible design 
space. If the weighting factors are small, the constraints exert little 
Influence on the move and the design may move in the direction of the con- 
straints . The behavior is dependent on the reduction factor applied to the 
weighting factors at each successive step, which is in turn related to the 
range of the weighting factors and the number of steps selected. A large 
number of steps (or a small reduction factor) can result in erratic behavior 
of the process due to the strongly repelling influence of the constraints . 

Too few steps, with the corresponding large reduction factors, may result in 
moves which impact one or more constraints before significant weight reduc- 
tions have been accomplished. 


The use of the approximate second derivatives in an adaptation of 
Newton's method results in a more efficient unconstrained minimization pro- 
cedure than does the use of the DFP algorithm (Reference 22). Since this 
latter procedure requires a nimiber of one-dimensional searches approximately 
equal to the number of design variables, and the number of such searches 
required with Newton's method is independent of the number of design vari- 
ables, the advantage of Newton's method increases as the number of design 
variables increases. As indicated in equation {6.kl) however, the matrix G 
must be inverted; this matrix is an n x n matrix where n is the number of 
design variables . Reference l6 states that this matrix is singular or 
ill-conditioned when the number of active constraints is small. To preclude 
the obvious difficulties which would otherwise result, a matrix G is used 
in place of the matrix G, the elements of which are shown in equation (6.U2). 


In this equation. 
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is the Kronecker delta and a value 


G. . = G. .(1 + e6. . ) 
ij ij ij 


(6.142) 


of G=0.01 is found to be satisfactory for the systems evaluated thus far. 


For each step in the optimization process, a one-dimensional search is 
conducted to determine the minimum of the modified objective function. The 
direction of the move vector is determined as previously described (equa- 
tion (6.4l)), and only the magnitude of the move vector is varied during 
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the search. The derivatives of the modified objective function need only be 
evaluated once during each step, but the modified objective function is 
evaluated for a number of values of the move vector magnitude sufficient to 
define the minimum of the function. These evaluations of the modified objec- 
tive function require flutter speed solutions, stress analyses and/or other 
procedures appropriate to the determination of the penalty terms . Although 
the process of evaluating the acceptability of a design relative to the con- 
straints is common to all optimization methods, the penalty function method, 
employing a one-dimensional search, requires at least three such evaluations 
per step, whereas methods which do not employ the one-dimensional search need 
only one such evaluation. These other methods, however, normally require a 
greater niomber of steps to achieve an acceptable level of optimization, and 
therefore usually require a larger number of flutter speed derivative deter- 
minations. Considering that determining the flutter speed derivatives 
requires determination of the flutter root and two characteristic vectors, 
it is estimated that in terms of computational operations (and thus cost) 
each step in the penalty function procedure, requiring the calculation of 
three flutter roots, is approximately equivalent to three steps of a pro- 
cedure such as the modified Simodynes method described in Section 6.3. 

This trade-off must be taken into account in any comparative evaluation of 
methods, and must be determined for a realistic design procedure. The 
principal advantage of the one-dimensional search is that it provides a 
specific criterion for the determination of step size, resulting (usually) 
in greater weight reductions per step than obtained with methods employing 
arbitrary step size. The determination of step size by the present procedure 
can readily be implemented as a fully automated computational subroutine. 

6.U.3 Assessment of the Method - Based on the description of the optimization 
method presented in Section 6. it. 1 and the discussion of the characteristics 
of the method in Section 6.4.2, it is concluded that the penalty function 
method (as that term is used here) is an efficient optimization procedure 
possessing characteristics not found in other methods. The unique means of 
handling design constraints is, of course, the most noteworthy of these. 

This feature, along with the use of the one-dimensional search to determine 
step size, results in an optimization procedure which is particularly well 
suited to use in a completely automated routine. Detailed specifications 
for such a procedure should be relatively easy to develop. The principal 
element required for the use of this method, over and above the requirements 
of other methods investigated, is the determination of the derivatives of the 
constraint quantities (other than flutter derivatives) with respect to the 
design variables . Since it should always be possible to evaluate the con- 
straint functions for specified values of the design variables, these deriva- 
tives may be obtained by finite difference techniques if no better method is 
available. In terms of computational efficiency, it is difficult to compare 
the penalty function method with the other four methods considered here. 

From the foregoing it is clear that the penalty function method does, in 
general, require more computations per step than do the Simodynes or Rudisill- 
Bhatla procedures. Whether or not the more efficient optimization step of the 
penalty function method overcomes that disadvantage can only be assessed in 
terms of a realistic design case. 



6.5 A Method of Feasible Directions 


The method of feasible directions is an approach to structural opti- 
mization using direct minimization of a constrained function. This is in 
contrast to penalty function methods, treated in Section 6.U, which convert 
the constrained design problem into a sequence of unconstrained minimizations 
of a modified objective function. 

The method discussed here is based primarily on the method of Gwin and 
McIntosh (Reference 17) > which is in tiirn a generalization of a method 
developed by Zoutendijk (Reference 2k). Addition background material is 
derived from Vanderplaats and Moses in Reference 25. It shoiild be noted that 
the discussion presented here is limited to those characteristics inherent 
to the feasible directions method itself; other particulars of the method 
presented in Reference 17 are not treated. Thus, such elements as the method 
of generating aerodynamic parameters, solution of the flutter equation and 
computation of flutter speed derivatives are considered to be separate from 
the method of feasible directions. These and other details of a complete 
optimization procedure are considered elsewhere in this report. 

6.5-1 Description of Method - The method of feasible directions generates a 
sequence of design changes, each of which is both feasible (does not violate 
active constraints) and usable (reduces total mass). Each resizing direction 
is followed until a new constraint is violated, an active constraint is 
re-encountered or the total mass is minimized. The process can be visualized 
with the help of Figure 6-1, which reproduces Figure 6 of Reference 17. In 
this figure, the minimum size constraints are represented by the horizontal 
and vertical boundaries, the minimum flutter speed constraint by the curved 
boundary and the constant weight contours by the diagonal straight lines . 

The starting point for the process is at point A, which is on the flutter 
speed constraint boundary; the design proceeds in a direction away from the 
flutter speed constraint and in a direction of decreasing weight until a con- 
straint boundary is encountered. At that point, a new direction is generated 
and a new move executed. This process is continued until a point B is reached 
which approximates an optimum design. Figure 6-2, also from Reference 17, 
illustrates the range of acceptable design directions starting from a point B 
on a general nonlinear constraint boundary. 

For a constraint equation of the form given in equation (6.43), the condi- 
tion that the direction is feasible is given by equation (6.44), where {Vh} 
is the gradient vector of the constraint with respect to the design variables, 

m^ , and {s} is the direction vector of the design variables. 


h < 0 


(6.43) 


[Sj {Vh} < 0 


(6.44) 
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In the present context, a feasible direction is one which does not 
violate an active constraint, as defined by a linear approximation of that 
constraint. The line segment BC in Figure 6-2 lies on the boundary of the 
feasible region. The condition that the direction is usable is given by 
equation ( 6 . U 5 ) 5 where {vw} is the gradient vector of the objective func- 
tion, i.e., the total structural weight. Any direction between lines BD and 
BC in Figure 6-2 is then in the usable-feasible sector. 

[sj {vw} <0 (6.1+5) 


The particular direction vector used in the- present iridthod is found such 
that a scalar /3 is maximized subject to the conditions expressed in equa- 
tion ( 6 .U 6 ). In this set of equations, B is an adjustment factor which 


(a) 

[Sj {Vh} + 

0/3 < 0 


(t) 

[sj {vw} + 

/3< 0 

(6.1+6) 

(c) 

The length 

of [_sj is bounded. 



controls the direction of j_sj within the usable-feasible region; large 
values of 6 force the direction away from the constraint and toward the 
usable boundary BD, while a value of 6=0 results in a move direction along 
line BC. For nonlinear constraints, such as a minimum flutter speed con- 
straint, an intermediate value of 6 is used which approximately bisects the 
usable feasible sector. Reference IT suggests a value of 0 = 1.0. It will 
be shown, however, that the effect of any particular value of 0 depends on 
the normalization of the constraint gradients and the weight gradient. For 
linear constraints, such as minimum size constraints, a value of 0 = 0 is 
used in order to produce a move direction parallel to the constraint boundary. 
The condition that the length of (_sj is bounded is usually accomplished by 
limiting the elements as shown in equation (6.^7). 


<1 i=l,2,....,n 


(6.1+7) 


Once the value of 0 is selected, the suboptimization problem indicated 
by equations (6.1+6) is transformed to a standard form and solved by means of 
the Simplex algorithm. Some of the details of this procedure are given in 
Reference 17, and a more comprehensive description may be found in Reference 
26 . The details of the Simplex algorithm will not be repeated here, but 
some aspects of the procedure are discussed in the following section. It 
should be noted that the Simplex algorithm was used to generate the direction 
vectors for the numerical evaluations presented in Appendix A of this report. 
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In applying the constraint conditions indicated by equation (6.46(a)), 
only those constraints which are considered to be active at a particular 
design point are included in the direction-finding process. A constraint is 
considered to be active if the design in question falls within a specified 
tolerance band t of the constraint boundary. This constraint tolerance 
band is arbitrarily specified and may vary with type of constraint. As 
stated earlier, any particilLar design step continues in the specified direc- 
tion iintil a constraint is violated; at that point, a correction step is 
taken back into the constraint tolerance band and a new direction generated. 

6 . 5.2 Principal Characteristics of the Method - The principal distinguishing 
feature of the method of feasible directions is, as the name implies, the 
direction generating procedure. Some aspects of this procedure, and of the 
method of establishing the step size, are discussed in the following. 

Some uncertainty exists as to the handling of the constraint tolerance 
band , e . Reference 17 indicates that an appropriate constraint tolerance 
is established for each constraint or constraint type, and that when a con- 
straint is violated, the step size is corrected such that the end-point of 
the step is midway in the tolerance band. Reference 25, however, recommends 
the use of a large tolerance band for the initial phases of the optimization 
so that the constraints will become active early in the optimization pro- 
cedure and will remain active. Combining these two approaches would not 
appear to result in an efficient resizing procedure, since a large constraint 
tolerance band would prevent a close approach to any constraint. For the 
idealized test case of Appendix A, no particular advantage can be discerned 
in maintaining a large "pad" on the constraints; this would almost certainly 
result in an increased number of resizing steps for the weight reductions 
shown in Section A. 6. It may be that such a procedure would be useful in 
the application of the method to more complex design problems as a means of 
avoiding convergence to localized minima, but it is felt that such usefulness 
would be rare in practical situations. Based on the idealized test case, it 
would appear that the constraint tolerance for minimum size constraints should 
be equal to zero. A minimum size constraint should not be active until the 
minimum size is reached, since otherwise no further reduction in that design 
variable would take place while the constraint remained active. For minimum 
flutter speed constraints, it would seem that the constraint tolerance band 
for defining an active constraint should indeed be fairly large. The flutter 
speed constraints then become active early in the design process, and are 
therefore effective in the efficient resizing of the design variables. The 
constraint tolerance band for determining the end-point of a particular design 
step, however, should be essentially equal to zero. 

Assuming a reasonable value of the adjustment factor, 0 , the next 

design step will be directed well into the feasible region, so that no useful 
purpose is served by originating the step any appreciable distance from the 
flutter constraint boundary. 

For other types of constraint, the treatment of the constraint tolerance 
band can be based on similar reasoning. 
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As indicated earlier, the reconmiended adjustment factor for a minimiun 
size constraint or other linear constraint is S = 0. This choice of adjust- 
ment factor results in subsequent design steps proceeding parallel to the 
constraint boundary if the design steps would otherwise result in continued 
reductions of the constrained design variable. For nonlinear constraints 
such as flutter speed constraints, a nonzero value of 9 must be used in 
order to force the design direction away from the constraint boundary. 
Otherwise, any finite move amplitude would violate the flutter speed con- 
straint. This characteristic of the parameter 9 results in the term 
"push-off" factor being applied to it in some references. The most efficient 
value of the "push-off" factor depends on the particular structural opti- 
mization under consideration, but a value of 9 = 1.0 is usually chosen. 
Reference 17 states that this value produces a vector direction which 
approximately bisects the usable-feasible sector. Examination of equations 
(6.h6) demonstrates, however, that the value 9 which accomplishes this is 
dependent on the normalization of {Vh} and {VW} . A larger value of the 
normalized constraint gradient corresponds to a smaller value of 0 , while 

a larger value of the normalized weight gradient corresponds to a larger 
value of 0 . In the niimerical evaluations described in Appendix A, the 
weight gradient was a unit column as a consequence of the choice of design 
variables, and the flutter speed constraint gradient was normalized such that 
the value of the average element is unity. A value of 0 = 1.0, in conjunc- 
tion with that normalization, produced approximately the desired result. The 
more usual procedure of normalizing the gradients on the largest element would 
have produced approximately twice as much "push-off". 


It is noted that Vanderplaats and Moses, Reference 25, recommend a 
variable "push-off" factor which is a function of the distance from the con- 
straint boundary. In the present notation, this is expressed in equation 


(6.ii8), where G is the constraint tolerance 
constraint function at the design step. 


and 


is 


the value of the j^^ 
is chosen as unity. 





(6.48) 


Since the constraint functions have negative values anywhere within the 
feasible region, the effect of this treatment of the "push-off" factors would 
be to drive the design to the constraint tolerance boundary, where the value 


of 


is equal to zero. Although this approach has the advantage of pro- 


viding a uniform treatment of the "push-off" factors, it would appear that 
the resultant excursions of the design in and out of the constraint tolerance 
band might well offset any advantage derived from this approach. 
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The direction vector for feasible directions procedures is usually 
obtained by the Simplex method. As noted in Reference 25 > "the resulting 
direction tends to point towards the corners of a hypercube in the design 

space (S? =1 or -l)." The direction vectors obtained in the numerical 

evaluations of Appendix A certainly exhibit the indicated tendencies . Refer- 
ence 25 goes on to propose an alternate direction vector formulation, based 
on imposing bounds on the total vector rather than the individual elements . 
Although this formulation has not been evaluated in detail, it is presumed 
that the resulting vector might be somewhat more efficient since it is 
subject to fewer constraints. 

Returning to the original formulation, some elementary considerations 
of equations (6.H6) and equation (6.kT) will indicate the reasons for the 
resulting direction vector characteristics, and will suggest a method of 
generating the direction vector without the use of the Simplex algorithm. 
Referring to equations (6.il6), assume for the moment that only the flutter 
constraint is active. It will be seen that maximum j3 can only occur when 
both equation (6.U6(a)) and equation (6.i|6(b)) are equal to zero. If now the 
value of 8 is taken as unity, equation ( 6 .U 9 ) results when j3 is maximum. 





( 6 .U 9 ) 


Returning to equations {6.k6) and observing that the elements of {Vh} 
are generally negative, it is not difficult to see that the addition of weight 
to the design variable having the greatest increase in flutter speed per 
pound, and removing weight from the design variable having the least increase 
in flutter speed per pound, tends to increase /8 . The restriction on the 

size of the elements expressed in equation ( 6 .U 7 ) limits these elements 

to +1 and -1 , respectively, while the objective of maximizing /3 insures 
that these limits will be reached. Continuing the same line of reasoning, it 
will be seen that the resiilting direction vector must have +1 elements for 
the design variables with the higher flutter speed derivatives, -1 elements 
for the design variables with the lower flutter speed derivatives, and one 
element of intermediate value in order to satisfy equation ( 6 .U 9 ). If one or 
more design variables in the negative group is limited by an active minimum 
sizing constraint, the corresponding elements in the direction vector are set 
to zero and equation (6.49) balanced as before. It should be noted that the 
ranking of the flutter derivatives must be on the basis of rate of change of 
flutter speed per pound of design variable, which is the ratio of the elements 
of {Vh} to the elements of {vw} . During the course of the numerical 
evaluations presented in Appendix A, it was found that a move vector direction 
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identical to that generated by the Simplex method could he derived hy the 

procedures described. When there are two active flutter constraints, V 

and , the values of -= + determine which elements of IS I are 

2 dm^ dm^ 

+1 and -1 ; two equations similar to ( 6 .U 9 ) will determine two elements 
of LsJ that have Intermediate values. This was also demonstrated, numer- 
ically, to lead to the same LsJ as the one generated by the Simplex method. 
It is expected this approach can be expanded to more active flutter con- 
straints . It is suggested that it is more straightforward and may be a more 
economical means of generating the direction vector. 

Several methods of step-size selection are proposed in Reference 17, 
the simplest of which is to continue to increase step size in the prescribed 
direction until a constraint violation occurs. At this point, a corrected 
distance to a point within the constraint tolerance band is found by linear 
interpolation of the appropriate constraint function values . 

As mentioned earlier, it should be relatively simple to determine the 
step size which would terminate exactly on the boundary of the nearest linear 
constraint, and in view of the fact that no further reduction of that design 
variable will take place while the constraint is active, no correction into 
the tolerance band would appear to be required. A more direct procedure might 
then be to determine the step-size to the nearest linear constraint and then 
check for violation of other constraints using that step-size. If constraint 
violations resulted, the step-size would be reduced until the critical con- 
straint was just satisfied. For minimiim flutter speed constraints, the use 
of some form of Incremented Flutter Analysis to solve for the step-size 
necessary to satisfy the flutter constraint should result in a substantial 
improvement in the step-size search procedure. 

6 . 5.3 Assessment of the Method - The method of feasible directions is similar 
to the penalty function methods in that the constraint functions influence 
both the direction of the design step and the step-size. Pull automation of 
the method is doc\imented in Reference 27. A wide range of constraint types 
can be accommodated, the only requirement being that it must be possible to 
evaluate both the constraint functions and constraint gradients for each 
design step. In concept, multiple flutter constraints can be included, 
although there may be practical difficulties to be overcome outside the 
optimization proper. 

The move direction vector resulting from the usual form of the equations 
appears to be rather crude, and it is felt that an alternate procedure, such 
as suggested in Reference 25, might result in a more efficient move vector. 

The treatment of constraint tolerances and determination of step-size seem 
overly complicated, and significant improvements in these areas should be 
possible . 


61 


Overall, the method seems to be quite representative of the direct 
methods of flutter optimization and competitive with other methods evaluated. 
The ntimerical evaluations reported in Appendix A show the method to be much 
better behaved than would be indicated by consideration of the move vector, 
and the rate of convergence demonstrated on the simplified test case is quite 
satisfactory. As a result, the method of feasible directions must be con- 
sidered a strong candidate for further evaluation in a realistic design 
environment . 


6.6 An Optimization Method Using Incremented Flutter Analysis 


Incremented Flutter Analysis was conceived as a method for efficiently 
determining design parameters for external stores satisfying a predetermined 
flutter speed requirement (Reference 4). Its capability of determining the 
value of general design variables such that the flutter speed has a given 
value makes it attractive as a tool in optimization with flutter constraints. 

A study was therefore initiated to assess the usefulness of Incremented 
Flutter Analysis as a tool in an optimization procedure. It was used in a 
heavily interactive approach to optimization for flutter, using the Computer 
Graphics system. The resulting method of optimization was demonstrated on a 
simulated design problem based on a subsonic transport wing (Appendix A, 
Section A. 7) and was also used on an actual design problem (arrow wing super- 
sonic transport) . 

During this development, the method of Incremented Flutter Analysis was 
generalized to be consistent with the needs in a complex optimization program 
(Reference l). 

A breadboard prototype of an automated computer program was developed 
and demonstrated on the same simulated design problem. 

In the following, the main features of the program and its present form 
are presented based on data in Reference l8. 

6.6.1 Main Features - The optimization method is a resizing routine "that 
minimizes total mass while maintaining the flutter speed exactly at a required 
value . 

Key features of the method are that the resizing column is allowed to 
change direction, without the need to recalculate flutter speed derivatives, 
during a one-dimensional minimization process in which the value of a scalar 
ot is determined that minimizes the total mass. 
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In the methods of optimization discussed in the preceding sections, a 
column of design variable increments (resizing column) |Am^ | is defined as; 




( 6 . 50 ) 


where defines a direction in design variable space and the scalar 

a magnitude. The resized design is related to the original design by 


{mfl} = {mjj+ (Amj) 


(6.51) 


In this method the resizing column is 


= |d^(Q!^)| 


(6.52) 


i.e., the direction of 


[AmJ) 1 


is a function of the scalar a. The scalar a 


is discussed in Section 6.6.2. During one resizing cycle, i.e., for one 
value of the superscript k, one set of derivatives of the flutter speed 




9V 

with respect to the design variables, , is used. The column matrix 

8V ^ 

is a function of as well as sizing constraints. During the 

i 

one-dimensional minimization of the total mass with Q! as a variable, the 
method of Incremented Flutter Analysis is used to maintain the flutter speed 
exactly at the desired value. 

6.6.2 Present Form of Program - The resizing column defined as 

the sum of a basic resizing column adjustment column 6 |a^| . 
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The elements of 


satisfy the equation 




( 6 . 53 ) 


and thus would correspond to a zero change in flutter speed if V were a 
linear function of the m^'s. The scalar 6 of the adjustment column 
is determined such that 




( 6 . 54 ) 


results exactly in a zero change in flutter speed. 

8V 

In the following the notation ~ — = V . is used.. 

dm. mi 
1 

The design variables are divided in an R group and a Q group, 
such that 



( 6 . 55 ) 


The division between the R and Q group lies within the positive range 

of V . . Thus all design variables for which V . < 0 are in the R group, 
mi mi 

The largest V . is identified by i = m and the smallest positive V . 

mi mi 

by i = sp. The division between the R and Q groups is defined by the 
inequality 

V /V . - 1 
^mm mi ^ 

mm' msp 

where VR is an empirical value. It was found, by numerical experimenta- 
tion, that VR = 0.3 is an acceptable value for the test case reported in 
Appendix A. For i = m the left hand side equals zero; for i = sp it 
equals 1. Thus the design variables for which the inequality (6.56) is 
satisfied belong to the R group. 

The algorithm in Reference l8 is based on removing mass from each design 
variable in the R group that is not at minimum size. Each mass removal is 
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individuELlly coupled with a change of all design variables in the Q group 
such that equation (6.53) is satisfied. Mass removal from a design variable 


in the R group for which > 0 requires addition of mass in the Q 
group. If < 0 mass removal in the R group is compensated by masi 
removal in the Q group in order to satisfy equation (6.53). 


The column matrix 




is the change in the design variables in 


the Q group due to removal of mass from design variable r in the R 
group. In Reference 18 the distribution 


C V 
q _ mq 

C V 
m mm 


(C > 0) 


m 


(6.57) 


is used if mass is to be added to the Q group (C >0 and V . „ > O) and 

m mien 


V 

m mq 


( 6 . 58 ) 


if mass is to be subtracted from the Q group (C <0 and V . „ < O). 

HI lElSn 

Equation (6.57) expresses that more is added to the design variables with 

the higher values of V .. Equation ( 6 . 58 ) expresses that more is subtracted 

mi 

from the design variables with the lower values of V . . 

mi 

The distribution for removal of mass in the R group used in Reference I 8 

is : 



( 6 . 59 ) 


The foregoing leads to a basic resizing column that is the sum of two 
columns that do not "overlap" and thus can be written as one column 


{‘=il = 



( 6 . 60 ) 
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The column matrix 



is given by equation (6.59); 



by: 


C 

q 


-c 



mr mq ’ 


mq 


-C 


eV >0 
mr 


V 
mr 

n V 
q mq 


^reV <0 

mr 


(6.61) 


where n^ is the number of design variables in the Q group. 

Equations (6.59) and (6.6l) define the basic resizing column 
the scalar C* defining a magnitude. 



with 


If the value of C* is such that equation (6.59) leads to the violation 
of minimum size constraints, that equation is only used for the elements that 
do not violate the size constraints. The other elements are given values 
corresponding to the minimum size constraints. 


Violation of a minimum size constraint by design variables in the Q 
group is expected to be infrequent. The program, however, has provisions to 
guard against such violation. 


The distribution of the adjustment column 6 j^ij discussed in 

Reference l8. In the numerical example in Appendix A, = 0 except at 

9V 

i = m (maximum is suggested, however, that = C. for icQ 

i . 

and A^ = 0 for ieR is expected to be a better choice. 


The value of the scalar 6 is determined by means of Incremented Flutter 
Analysis. The determinantal flutter equation (according to equation (3.^)) is 
written as : 


d {(7 + i )k,g,V,p,m^+C^,6A^) = 0 (6.62) 

In equation (6.62), 7=0, g=0, V and p have given values; m^ 

are the values of the design variables at the beginning of the current 

resizing step; satisfies equation (6.53). Equation (6.62) is solved by 

two-dimensional Regula Falsi for k and 6 . The complete resizing column is 
then determined by equation (6.5^). 
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If no size limitations become active, simple substitution of equa- 
tion (6,59) into equation (6.6l) shows that all elements of the resizing 

column proportional to C*. Different values of C* can be 

C^l + computed and M = |^lj|m^ + Arn^j can be 

computed as a function of C*. This is, in principle, the one dimensional 
minimization that determines the direction as well as the magnitude of the 
resizing column. The total mass M does have a minimum due to nonlinear 
effects, which are taken into account by means of Incremented Flutter Analysis 

When size limitations are active not all elements of |c^ 

tional to C* and nonalgebraic operations are needed to obtain M as a 
function of C*. The logic for this is presented in Reference l8. 

To facilitate initiation of the one dimensional minimization, a variable 

a = - Z C is used. It is the total mass removed from design variables in 
r r 

the R group and has a very simple relation to C* . The quantity a has a 
simple physical meaning, independent of the number of design variables. An 
initial value can be chosen as a fraction of the total mass represented by 
the design variables. 

The numerical example of Reference l8 is presented in the Appendix as 
Table A-11. It suggests that the method is very powerful in reducing the 
total mass by a large fraction (>80^) of the difference between the current 
total mass and the minimum mass in a single step. 

Although it is not essential to the overall method, it should be noted 
that the program as described in Reference l8 uses the method of Incremented 
Flutter Analysis to determine the values of 9V/9m^ by means of a finite 

difference approach (Reference l). 

6.6.3 Concluding Remarks - The approach taken in this method is distinctly 
different from each of the other methods discussed, although it contains 
elements of several of these methods. One distinguishing feature is that the 
flutter speed is held exactly at the constraint value. In fact, this feature 
is used to include the nonlinear character of the flutter speed as a function 
of the design variables and makes it possible to do a one-dimensional minimi- 
zation of the objective function itself, rather than of a modified objective 
function as in the penalty function method. 

Another distinguishing feature is the departure from the traditional dis- 
tributions in the resizing column, which are largely based on the gradient of 
the velocity and the total mass. It has been demonstrated that empirically 
generated distributions can lead to rapidly converging optimization procedures 


are propor- 


assumed 





A third feature is the use of Incremented ' Flutter Analysis. It is not 
considered advantageous to use Incremented Flutter Analysis to determine the 
derivatives of the flutter speed by the finite difference method. However, 
Incremented Flutter Analysis, as executed with the help of the two-dimensional 
Regula Falsi method for solving two non l inear equations with two unknowns, 
is an efficient method for keeping the flutter speed exactly constant. This 
use of Incremented Flutter Analysis could be used advantageously in some of 
the other methods of optimization: keeping the flutter speed exactly constant 
facilitates the observance of a convergence criterion for minimum total mass 
since side effects due to drift of the flutter speed are avoided. 


6.7 Comparison of Optimization Methods 


6 . 7.1 General - In comparing the optimization methods discussed in 
Sections 6.2 - 6 . 6 , many differences in procedural detail are apparent. 
Specifically, differences in generating the distribution and magnitude of 
the resizing column can be recognized. Of these two, however, the more basic 
difference relates to the determination of the magnitude of the resizing 
column, or step size. Several of the methods - the weight gradient optimiza- 
tion of Simodynes (Section 6.3) and the velocity gradient, mass gradient, and 
gradient projection methods of Rudisill-Bhatia (Section 6 . 2 ) - employ arbitrary 
step sizes. In contrast, the penalty function method (Section 6 . 4 ), the 
method of feasible directions (Section 6.5) and the method incorporating 
Incremented Flutter Analysis (Section 6.6) all make use of a step size deter- 
mined by well defined criteria. In the following, these two groups are dis- 
cussed separately and then a candidate resizing procedure is synthesized, 
based on the analytical and numerical evaluations conducted thus far. 

6.7.2 Arbitrary Step-Size Procedures - The resizing procedures employing an 
arbitrary step size are characterized by a well-defined resizing cycle which 
is usually simple and straightforward. In general, one flutter solution 
(with two characteristic vectors) and one set of flutter derivatives are 
required for each step. The procedures based on constant flutter speed 
(Simodyne's weight gradient and the gradient projection search of Rudisill- 
Bhatia) rely on a linearization, based on the flutter velocity derivatives, 
to hold the flutter speed constant. As a result, the actual flutter speed 
tends to drift (downward, in most practical resizing exercises) and must 
periodically be corrected. It is this tendency which effectively limits 
step-size, since large excursions from the required flutter speed, are unde- 
sirable. Rather than attempting to maximize step-size, however, a moderate 
step-size is chosen, based on experience and engineering Judgment, and the 
attendant penalty of an Increased number of steps is accepted. 

In terms of specific procedures, the Rudisill-Bhatia methods should be 
somewhat more efficient than the weight gradient method of Simodynes . As 
indicated in Section 6 . 3 , the frequency constraint imposed by this latter 
procedure results in a degree of approximation which probably cannot be 



Justified on the basis of the resulting simplifications. In the modified 
procediire (Section 6.3.3) > this frequency constraint is removed and the 
resulting procedure is shown to be similar to the gradient projection 
search of Rudisill-Bhatia. In comparing these two procedures, the Rudisill- 
Bhatia approach has a significant advantage in that it does not require the 
selection of a dependent design variable, and thus eliminates the resulting 
influence of this choice on the performance of the procedure. The use of 
the flutter velocity derivatives as developed by Rudisill-Bhatia, rather than 
the normalized derivatives of Simodynes, has some advantage in a procedure 
which incorporates a nonzero flutter velocity increment. This would be the 
case when using the flutter velocity gradient to define a design variable 
distribution to increase the flutter speed of an initially deficient system, 
or to make small flutter velocity corrections during the resizing cycles. 

For use in constant flutter velocity procedures, however, it is shown in 
Section A. 3. 3, Appendix A, that the two forms of the derivative may be used 
interchangeably . 

The two most useful procedures employing arbitrary step size would then 
appear to be the velocity gradient search and the gradient projection search 
of Rudisill-Bhatia. As envisioned in Section 6 . 7 .^, this former procedure 
would not be implemented using an arbitrary step-size. The most useful 
application of the velocity gradient search procedure is in increasing the 
flutter speed of a flutter-deficient system to the required flutter speed in 
a nearly optimum manner. It is shown in Section A.U.l, however, that the use 
of a single step to accomplish this is probably the most effective procedure. 
As a consequence, it seems reasonable to use the Incremented Flutter Analysis 
technique (cf. page 43) to determine the step-size necessary to satisfy the 
flutter speed constraint exactly. The gradient projection search could also 
be improved by a modification of the resizing column such that the mass 
gradient component is replaced by the reciprocal of the flutter speed 
derivatives, resulting in a resizing column defined by equations (A.3)» (A. 4) 
and (a. 5 ) of Appendix A. A comparison of Tables A-5 and A-8 of Appendix A 
indicates that, at least for the idealized test case evaluated there, a 
significant Increase in efficiency results from this modification. 

6 . 7.3 Defined Step-Size Procedures - In resizing procedures employing a 
defined step-size, an attempt is made to maximize the step-size so as to 
derive the maximm benefit from a single resizing step. This maximum step- 
size is determined by a well-defined set of criteria, usually involving the 
condition of the current design with respect to the design constraints. 
Evaluating these criteria usually involves the determination of the flutter 
speed, along with other constraint conditions, at several points along the 
move path during one step. In contrast to the arbitrary step-size procedures, 
then, the defined step-size procedures normally result in a greater mass 
reduction for each set of flutter derivatives calculated, but at the expense 
of a greater number of required flutter solutions. 

Each of the three defined step-size procedures discussed in the preced- 
ing sections has distinct characteristics, and meaningful comparisons are 
difficult to make. Some general observations are possible, however, allowing 
some tentative conclusions to be drawn. 
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The penalty function procedure (Section 6 .k) is perhaps the most versatile 
of the three procedures considered. The treatment of constraints is straight- 
forward, making automation of the procedure particiilarly simple. A step is . 
terminated before a constraint violation occurs, but the constraints are 
continuously active and exert some influence on the direction of each resizing 
step. The efficiency of the method is to a significant degree dependent on 
the handling of the penalty term weighting factors associated with the con- 
straints (equation (6.38)), so that some Judgement and experience are both 
required in the selection of these factors. A more troublesome difficulty 
might be encountered in the use of the direction generating algorithm based 
on Newton's method. As shown in equation (6.I1), the matrix of coefficients 
of the design variable second derivatives must be inverted, and it is 
indicated in Reference 16 that this matrix may be singular or .ill-conditioned. 
Although means to avoid this problem are suggested, computational difficiilties 
may still arise in practical design tasks involving large numbers of design 
variables . 

The method of feasible directions (Section 6.5) provides a treatment of 
constraints which is different from that of either of the other two methods 
considered here. For the flutter speed constraint (and presumably other 
nonlinear constraints) the approach appears to be quite satisfactory. It is 
similar to that of the penalty function method, in that the "push-off" factor 
can be considered as analogous to the penalty weighting factors of that 
procedure. In the feasible directions method, however, the move direction is 
restricted to the usable as well as feasible region. As a consequence, each 
move must reduce the objective function (total mass) as well as avoid a 
violation of the design constraints. In contrast, the penalty function move 
must reduce the modified objective function, but not necessarily the objective 
function itself. For linear constraints, such as minimum sizing constraints, 
the feasible directions method is formulated such that the constraints are not 
active until a constraint violation occurs, and the resizing step is 
terminated at that point. Prior to the constraint violation, the constraint 
exerts no influence on the move direction and therefore such constraint 
violations are normal occurrences. For the idealized test case evaluated 
in Appendix A, this characteristic did not appreciably degrade the efficiency 
of the procedure; only four sizing constraints became active dviring the 
optimization, and a significant weight reduction resulted from each step. In 
a more realistic design case, with a large niomber of minimum size constraints, 
it is anticipated that a significant number of short, ineffective moves would 
result from sizing constraint encounters. Once a linear constraint becomes 
active, however, the constraint is incorporated in the move direction so 
that subsequent moves take place along the constraint boundary. The conditions 
imposed on the direction of the move vector result in design variable incre- 
ments that, in general, are equal in magnitude, are positive for the design 
variables with the higher values of 9 V/ 9 m and are negative for the design 
variables with the lower values of 9 V/ 9 m. Thus there appears a lack of 
differentiation between the design variable increments with equal algebraic 
sign. It should be noted, however, that the results of the numerical evalua- 
tions described in Appendix A do not substantiate this lack of efficiency. 
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The procedure incorporating the use of Incremented Flutter Analysis in 
a formalized resizing procedure (Section 6.6) utilizes a concept to define 
step-size which is significantly different from that of either of the other 
two methods discussed here. In principle, a distribution of design variable 
increments is defined which reduces total mass and which, on a linear basis, 
produces no change in flutter speed. Since flutter speed is not a linear • 
function of the design variables, however, any finite value of this incremental 
distribution will produce (in all practical cases) a decrease in flutter 
speed. For several values of the incremental distribution, the value of an 
adjustment increment is determined by the Incremented Flutter Analysis tech- 
nique which brings the flutter speed exactly back to the required value. For 
some value of the incremental distribution the total mass will be a minimum; 
this point defines' the end of the step. By this means,. the nonlinearities in 
the flutter constraint are explicitly accounted for, and the maximum mass 
reduction for a given set of values of the flutter speed derivatives is 
achieved. In the process of determining the step-size associated with the 
minlmiira mass, minimiim size constraints are enforced as necessary. To that 
extent, the direction vector of the design variables is a function of the 
move amplitude, the direction conforming to the size constraints. 

As presented in Reference l8 and as used in the numerical evaluations 
of Appendix A, the method described in Section 6.6 differs in two other 
respects from the other methods evaluated. The flutter derivatives are 
obtained through the use of Incremented Flutter Analysis in the form of 
increments in each individual design variable required for a reference change 
in flutter speed. The results of the numerical evaluation indicate that this 
finite difference form of the flutter derivatives results in values compa- 
rable to those obtained from the analytic form, and that the two forms may be 
used interchangeably. It is recognized, however, that the use of this pro- 
cedure for obtaining the flutter derivatives - which requires the equivalent 
of a flutter solution for each design variable - is not efficient for a 
practical design task involving a large number of design variables. In such 
a case, the use of the analytic form of the derivative would be more econom- 
ical. The other area in which this method differs from those previously 
discussed is in the formation of the resizing move vector. The separation 
of the design variables into two groups, those with higher values of 9V/9m 
and those with lower values of 9V/9m, is empirical, and it is not clear 
that the criterion used in Section 6.6 would be efficient in all cases. The 
distribution of the increments for the design variables with the higher values' 
of 9V/9m is proportional to the velocity gradient, but the distribution for 
the design variables with the lower values of 9V/9m is a second order 
function of the reciprocals of the velocity derivatives. The results of a 
numerical evaluation using the move vector of Section A. 3. 3, Appendix A, in 
place of the move vector described in Section 6.6 indicate that the effi- 
ciency of the two move vectors is approximately equal, at least in terms of 
the idealized test case of Appendix A. In view of this, it is considered 
that the more complex move vector presented in Section 6.6 is not justified 
on the basis of present results. 
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6,'J.k Formulation of a Resizing Procedure - Based on the evaluations 
presented in this section and the results of the numerical evaluations 
presented in Appendix A, a resizing procedure can he formulated which will 
result in an improved performance over that of any of the specific methods 
discussed."' It is recognized that the evaluations of the resizing procedures 
are not complete; any promising procedure must he further evaluated in terms 
of a realistic design task in order to arrive at firm conclusions. In 
particular, further evidence must he obtained to determine the relative 
efficiency of the arbitrary step-size and defined step-size procedures. 

For the time being, however, it will he premised that the problems associated 
with a practical design task will dictate the use of a defined step-size 
procedure. These problems, some of which are discussed in Section 7, would 
seem to indicate that the greatest possible mass reduction shoTild be obtained 
for each step, since the structural reanalysis required per resizing step may 
be much more extensive than is generally recognized. 

As qualified by the preceding paragraph, the preferred resizing pro- 
cedure may be described in terms of the following characteristics: 

Flutter speed derivatives are of the analytic type, calculated by the 
method of Rudlsill-Bhatia (Section 6.2), possibly generalized by taking into 
account the derivatives of the vibration modes with respect to the design 
variables . 

The initial resizing step is one which increases the flutter speed of 
the flutter-deficient design to the required value. This is done in a 
nearly-optimum manner by the addition of design variable increments distributed 
according to the velocity gradient. The total required flutter speed incre- 
ment is obtained in a single step, and a form of Incremented Flutter Analysis 
is used to determine the magnitude of the step required to satisfy the flutter 
constraint exactly. 

Subsequent resizing steps are performed at constant flutter speed, using 
the technique of minimization of the objective function (mass) described in 
Section 6.6 in order to determine the step-size. As discussed in that section, 
the primary distribution of design variable increments is such as to produce 
zero flutter velocity change on a linearized basis. Incremented Flutter 
Analysis is then used to determine the magnitude of an adjustment colvann of 
increments required to maintain the actual flutter speed constant. The 
total mass of the design variables, including both the primary and adjustment 
distributions, is determined as a function of step-size, and the step-size 
corresponding to minimum mass chosen. Using this configuration as a starting 
point, the resizing cycle is repeated. 

The sizing constraints are satisfied in the manner of Section 6.6, with 
the direction of the move vector being modified as constraints are encountered 
during the minimization of the objective function. Note that the flutter 
constraint is satisfied at each substep of the minimization. 

The move vector for design variables corresponding to positive values 
of BV/9m is based on the combination of the velocity gradient and the 
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negative reciprocals of the flutter derivatives shown in Section A. 3. 3. 

For design variables corresponding to negative values of 0V/9m, modified 
distributions will be used, some options of which are discussed in Reference 1 

, t 

A s\ammary discussion of this procedure and an example of numerical 
results are presented in Reference 28. 


7. CONSIDERATIONS RELATED TO A REALISTIC 
DESIGN ENVIRONMENT 

In the literature on optimization with flutter constraints , the methods 
presented are illustrated with examples of varying complexity. One example 
in Reference l6 is based on 156 structural degrees of freedom and 23 design 
variables. The example in Reference 29 is based on 150 degrees of freedom 
and 100 design variables. As much as these numbers surpass the corresponding 
numbers in earlier examples in the literature, they fall short of what may be 
encountered in a realistic design environment. Thus, problems that may result 
from such an environment remain unexposed. During the present work, several 
aspects of dealing with an actual design have been examined. They are dis- 
cussed in the following sections, together with other aspects to which little 
or no attention could be given. 


7.1 Structural Model 

The mathematical model representing the structure obviously is an 
important element of structural optimization with flutter constraints. Finite 
element structural models with thousands of elements and a corresponding 
number of nodal displacements as degrees of freedom are used for stress and 
stiffness analysis of a given structure. A duplication of effort can be 
avoided if the same structiiral model can be used for flutter optimization. 

It should be noted that for a flutter analysis, or a loads analysis 
including aeroelastic effects, the refinement of a multi-thousand element 
structural model is not required. For the ciirrent type of subsonic transports 
a. relatively simple beam model suffices for flutter. For supersonic trans- 
ports, however, as are in existence and projected for the future, a simple 
beam model is inadequate and a finite element model must be used. This 
implies that methods of optimization for flutter must be able to handle finite 
element type structural representation. 

The typical structural model for stress analysis has a number of degrees 
of freedom that exceeds what at present seems practicable for the repetitive 
vibration analysis expected in a flutter optimization program. A reduced 
ntanber of degrees of freedom can be obtained by using a stiffness matrix of 
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reduced size in the vibration analysis or by generating a coarser finite 
element model for aeroelastic analyses, which may or may not require further 
size reduction. 


The common approach to coordinate reduction is one in which coordinates 
to be eliminated are assimed to have zero loads. This is often called static 
reduction in contrast with the approach of Reference 30; which can be called 
dynamic reduction, since the reduction is a function of the frequency. If 


the basic stiffness matrix 
the reduced matrix, * 



is partitioned as indicated in equation (7-l), 


is given by equation (7.2). 



(7.1) 


H = [■'ll] - 


If the incremental stiffness due to increasing the design variable 
a unit amount is |aK^J, the basic stiffness matrix, as a function of the 
/3^'s is: 

[Vi’i’] = [V“>] K] '■f 


where is defined relative to a reference value. 

In general, therefore, ^ nonlinear function of due to the 

triple product and inversion in equation (7.2). Thus, to compute 
the coordinate reduction represented by equation (7.2) must be repeated for 
each combination of values of B. . 

In most flutter optimization procedures, the derivative of the stiffness 


matrix, respect to many design variables B^ is required. 

The procedure is as follows: 

Equation (7.2) is equivalent to 


H H H 


( 7 .^) 
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where 



(7.5) 


The derivative of the reduced stiffness matrix with respect to any 
design variable ( 3 ^, evaluated at a given combination of values of is 

then defined by: 

[akJ[0E(?,)] (T.6) 


where AK . 

X 


is consistent with equation (7.3). 


The application of equation (7-6) is as follows. The incremental stiff- 


ness matrices 


N 


are invariant during the optimization process. As a 


new set of design variables, /3^, is defined during a step in the optimiza- 
tion process, is computed according to equation (7*3) and accord- 
ing to equation (7.2). The reduced stiffness matrix can then be used 

in a vibration analysis . The associated coordinate reduction matrix [gr] 
is formed and used in the triple product of equation (7.6) to compute the 


derivatives of 


[kJ Vi 


th respect to all design variables. 


When the number of structural coordinates is not too high, the coordinates 
to be eliminated can be restricted to those that are of no interest to aero- 
elastic analyses and, in fact, can be considered unloaded. However, when the 
structural model is designed for stress analysis it may be desirable to 
eliminate coordinates that are of interest to aeroelastic analyses, such as 
deflections perpendicular to lifting surfaces. This usually means that 
coordinates must be eliminated that have associated inertia. If that is the 
case, equations (7«^) and (7.6) must be applied to the mass matrix as well 
(Reference 31). 


In using a coarse grid finite element model for aeroelastic analysis, the 
aim is to reduce the number of coordinates to be eliminated to a minimum. In 
a finite difference approach, as described in Reference 32, the only struc- 
tiiral degrees of freedom are deflections perpendicular to the lifting siirface. 
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Which approach will he favored in future optimization work is hard to 
foresee. There seem to he three areas of investigation that could lead to 
significant development. 

It seems most advisable, because of the directness of the approach, to 

speed up the computation of >] for arbitrary sets of in the basic 

finite element analysis system. Possibly approximate methods can be developed, 
which are valid for a few resizing steps,, after which an exact updating takes 
place. It seems self-evident that the last updating in an optimization should 
be exact. 


A second area of investigation could be based on an approach used with 
some success at the Lockheed-California Company. In it the reduced stiff- 
ness matrix is approximated by a polynomial function of the design variables: 




(7.7) 


where the summation is over i = 1 -* n and j = !-♦ n. 


Such a polynomial can be an acceptable approximation over limited ranges 
of the values of the design variables. Since the stiffness is represented as 
an explicit function of the design variables, it can readily be evaluated for 
any arbitrary combination of values of /3. . The derivative of the stiffness 
matrix is : 


One element of 


air = M ^ 


(7.8) 


is approximated by 


K (/3.) SK (0) + S/3.A. + Z0./3.B., 

r 1 r 11 '^11 ij 


(7.9) 


To determine the values of the coefficients A. and B_, K ( /3. ) 

1 ij r 1 

2 

must be computed for n+n values of /J^. Thus, to define the polynomial 
expression in equation (7.7) it is necessary to compute 

1+n+n^ linearly independent columns |^i| ’ the help of equations (7.3) 

and (7.2) . 
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On an arrow wing, where design variables were torsional and bending 
stiffness over certain areas of the outer wing, it was found there was little 
coupling between the design variables. When that is the case equation (j.j) 
can be written as : 

= [k,(0>] ^ (7.10) 

Only l+2n evaluations of are necessaiy.to compute ( 0 )] 

and all the coefficient matrices |aJ and equation (T.IO). 

A third area of investigation is the use of the aeroelastic model. In 
that case it seems mandatory that a direct two-way relationship be developed 
between the sizing in the stress model and the sizing in the aeroelastic 
model . 


7.2 Multiple Flutter Speed Constraints 


Although the problem of multiple flutter speed constraints is addressed 
in the literature (e.g.. Reference 33 )» examples in the literature, used to 
illustrate methods of optimization for flutter, are all restricted to one 
flutter speed. Often it is indeed possible to eliminate all flutter con- 
straint violations by eliminating the most critical flutter speed. In gen- 
eral, however, the possibility of more than one active flutter speed constraint 
must be anticipated. 

Formally, the penalty function method (Section 6.U) and the method of 
feasible directions (Section 6.5) have built-in capability to handle multiple 
flutter speed constraints. The other methods discussed in Section 6 require 
added logic to handle multiple flutter speed constraints. 

Reference 3^ makes use of the multi-constraint capability of the penalty 
function method by requiring that the flutter roots for selected values of 
the reduced frequency, k, correspond to combinations of speed and damping 
that provide adequate damping within the flight envelope (see also Section 7-3). 


The optimality criterion for one flutter speed constraint is : 

3V ^ 8V 

9m. 9m, 

1 J 


(7.11) 


where i and j refer to free design variables, i.e., design variables that 
are not at a sizing constraint (References 29 and l). 
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For two flutter speed constraints the optimality criterion is 
(Reference 1): 


1 

1 

1 

CD 


f!i 

am^ 

9^2 

9m^ 


9^ 

9V2 

am^ 

9m^ 

ams 


( 7 . 12 ) 


Equation (7.12) must be satisfied for any combination of three free 
design variables. It is satisfied if: 


8v^ aVg av^ av^ 

am. am. am. am. 

1 1 J J 


(7.13) 


Extension of this criterion to more than two flutter speed constraints 
is straightforward. 

av^ avg 

It was found that the value — determines the resizing column 

i i 

that is generated by the feasible direction method of Reference 17 (Sec- 
tion 6.5) with a push-off factor 0=1. It is believed that this value can 
also be used in defining a resizing col\amn if the optimization is based on 
the methods discussed in Sections 6.2, 6.3 and 6.6. This is further discussed 
in Reference 1. 

To demonstrate multiple flutter speed constraint capability, the nximerical 
examples must relate to a realistic design environment in which two or more 
in-flight modes lead to flutter speeds below the minimum required flutter 
speed. These in-flight modes may be unrelated flutter modes for a particular 
weight configuration of the airplane and one particular Mach number, or they 
may be related or unrelated flutter modes for more than one weight configura- 
tion and Mach number. 


7 . 3 Damping Constraints 

In the discussions in Section 6, the emphasis is on flutter speed con- 
straints . This is in recognition of the flutter speed margins as defined by 
the Federal Airworthiness Regulations and military specifications. The 
requirements of the Federal Airworthiness Regulations are illustrated in 


78 



Figure 7-1* The airplane shall be designed to be flutter free within the 

altitude-speed envelope defined by 1.2 and h = -3100 m 

(-10,200 ft). In addition to this flutter speed requirement there is an 
implied requirement for adequate modal damping within this envelope. This 
is illustrated in Figure 1 - 2 : for all flight conditions id. thin the flight 

envelope there must be a certain amount of positive damping. In k-method 

terminology this means g < in p-k-method terminology 7 < ' 

Both g and 7 are negative quantities. From the definitions of g 

1U3>X 1113^ 

and 7 it follows that for small values g s 2^. 

To satisfy the general damping constraint, the inequality constraint 

Y <7 (or g<g ) must be invoked at several speeds below 1.2 V„, 

' ~ 'max “ max ^ D 

for all in-flight modes of interest, for several Mach numbers and for the 
airplane weight configurations to be considered. In a realistic design 
environment this may lead to hundreds of inequality constraints . It is 
obvious there are practical difficulties associated with that many constraints. 

Experience shows that usually only very few damping constraints are 
active and they are associated with hump modes. Sometimes damping constraint 
violations by hump modes disappear as the structure is resized to eliminate 
the most critical flutter speed violation( s ) . In anticipation of the need 
for invoking a minimum hirnip mode damping constraint. Section 3-^ presents a 
procedure to directly determine the minimum damping of a hump mode. 

The optimality criterion for an active damping constraint is : 


97 = 97 

9m. 9m. 

1 J 


(T.lit) 


where i and j refer to free design variables. It corresponds to the 
optimality criterion for a flutter speed constraint, the derivation of which 
(given in Reference l) can be generalized to arbitrary constraints. 


For a combined flutter speed and damping constraint, the optimality 
criterion is : 


111 

9V 9V _9V 

9m^ 9m^ 9m^ 

97 97 97 

7m^ rm^ 7m3 


(T.15) 
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Equation (7.15) must be satisfied for any combination of three free 
design variables. It is satisfied if 


^ + r-§21 =-^ , r 97 
dm. dm. dm dm. 

1 1 J J 


(T.i6) 


Here C is an arbitrary constant that can be used to create compatible 
units or to assign a different weighting to the two constraints . 


In several methods of optimization this criterion can provide a guide 
towards generating a resizing column. 

Further development of a practical method. of including damping con- 
straints in flutter optimization should be based on mmierical examples in a 
realistic design environment. 


7.^ Mass Ballast 


The literature pays little or no explicit attention to ballast (dead 
weight) as a design variable. The reason for this omission is understandable: 
any method of optimization that can handle design variables representing 
related stiffness and mass changes can handle a design variable representing 
a mass change only. 

Adding mass ballast may be a more efficient way of raising the flutter 
speed to its required value than structxaral stiffening. That appeared to be 
the case on one of the United States supersonic transport designs and in the 
example treated in Reference 29. In the latter case, however, it is not clear 
whether the modified strength requirements due to the addition of ballast are 
accounted for. Reference 29 demonstrates a potential complication associated 
with mass ballast as a design variable: as ballast is added in a particular 

region, the flutter speed derivative with respect to the mass ballast changes 
from negative to positive. As stated in Reference 29 , if this phenomenon 
occurs, an automated resizing procedure may fail to recognize the beneficial 
effect of a larger amount of ballast since an infinitesimal amount of ballast 
proved to lower the flutter speed. Until this aspect of flutter optimization 
has received more attention, considerable engineering judgment should be used 
in handling mass ballast as design variables. 

In view of the preceding paragraph, it may be convenient to first consider 
stiffness design variables and their associated masses only for raising the 
most critical flutter speed to the desired value and the subsequent optimiza- 
tion at constant flutter speed. Starting with the original deficient con- 
figuration, it is then determined whether any mass change without stiffness 
change is more efficient than the most efficient stiffness change in raising 
the flutter speed to the desired value. Incremented Flutter Analysis can be 
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used to directly determine the amount of ballast needed to meet the flutter 
constraint, regardless of any changes in sign of the flutter speed derivatives 
as a function of the amount of mass ballast. If mass ballast is more effec- 
tive than optimum stiffening, the mass ballast design variables of interest 
can be added as design variables for a final optimization process. 


7.5 Interface With Strength Optimization 


Combined optimization for flutter and stress has been demonstrated with 
simple structural models and/or under simplifying assxamptions (References 9, 
16 and 29) . 

In References 9 and 16 the penalty function method is used and, in order 
to reduce the number of stress constraint derivatives to be evaluated, the 
stress constraint is reduced to one constraint per loading condition. The 
effect of this can be easily seen in the following formulation of the penalty 
function approach. 


The modified objective function may be represented by: 
n n n 


0(m. ) = 2 m. + 


i=l 


r r r 

'"'“i' - ''r j=i 1=1 "1 - “i 


(7.17) 


where: 


V(m^) 


O’, (m. ) 

J 1 


m. 

1 


m. 

1 


flutter speed 

minimum desired flutter speed 

stress in jth element j =1 ... n. 

maximuim allowable stress in jth element 

design variable; mass associated with ith design element; 

i = 1 ... n. 


minimum allowable value of m. 

1 




m 


penalty weighting factors 
It is understood that in equation (7-17) there are as many 


V(„.) - Vj, 


terms as there are flutter speed constraints, and as many ^ •= — — 7 ^ 

._T cr. = O', (m.) 


. 1 * 

J=1 J 


J 1 
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terms as there are design load conditions. For this discussion it is suffi- 
cient to assiame one flutter speed constraint and one design load condition. 

As a part of the determination of the resizing column, the partial 

derivatives 90(m. )/9m. are used. 

1 1 


90 

9m. 


= 1 - 


V 


(v(mi) - V^) 


n 

2 9m. 9m. . - a. - cr. (m.) 

1 X J=1 J .11 


m 


("i -”i)^ 


(7.18) 


9V 

For each design variable, one derivative — and n derivatives 

r ^“i 

9 """o* 

S — ■= 7 T must be evaluated. Thus a total of n flutter speed 

9m. cr. - cr.Cm.) ^ 

1 J J 1 


derivatives and n stress derivatives must be evaluated. In References 9 


and l6 the mmiber of stress derivatives is reduced to n by evaluating the 
derivatives by means of finite differences but performing the differentia- 
tion after the summation; i.e., the following identity 



n 


0-^ - o-j(-i) 


= s 


J=1 (tj - 


9d. 

0_ 

2 9m^ 


(7.19) 


which requires the evaluation of n derivatives per single design variable, 
is replaced by the finite difference representation: 


n 

A. S 

'‘"i 0=1 ^0 - ''o<"‘i> 


6m. 

1 


n 

U=1 "j 


^cr 


n 


cr. (m. + 6m. ) 
J 1 1 


- 2 ^ 


‘0- 


•_n cr. - (r.(m. ) 
j=i 0 a 1 


( 7 . 20 ) 


requiring the evaluation of only one derivative per single design variable. 

This substitution does not affect the one-dimensional minimization that is 
part of the method used in Reference l4, but it does affect the direction 
of the resizing vector. Instead of each elemental stress contributing 

90 

individually to one contribution representing an average stress 

i 

penalty term is used. No studies investigating the effect of this substitution 
have been reported. 
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In Reference 28, the optimality criterion = constant for all i is 

i 

used in the flutter optimization and the fully-stressed-design criterion for 
strength optimization. The two optimizations are perfonned alternately until 
a converged design is obtained. 

The fully-stressed-design criterion does not, in general, lead to a 
minimxim weight structure, at least not for a redundant structure (Refer- 
ence 35). Furthermore, without further investigation there is little ground 
for expecting that alternating flutter and strength optimization will lead 
to a converged design if both flutter constraints and stress constraints 
are active. This leads to the conclusion that ideally flutter and strength 
optimization should take place simiiltaneously as is done in References 9 
and l6. 

It woTild seem that the adequacy of the methods of References 9 » 
and 29 has not been demonstrated when applied to a practical design problem. 
Conceptually the methods of References 9 and l6 allow as many stress con- 
straints per load condition as there are independent stress constraints 
(which may be larger than the number of elements ) . Thus the n\miber of inde- 
pendent stress constraints can be very large. This has led to the definition 
of one stress constraint per load condition, which, however, removes the 
possibility of independent stress constraints contributing individually to 
the move vector direction. Possibly other composite stress constraints can 
be defined such that the mmiber of derivatives to be evaluated is reduced 
while retaining the contribution of each constraint to the move vector 
direction. 

Additional investigations in the areas of structural modeling and combined 
optimization for flutter and strength are needed before conclusions regarding 
the best approach can be formulated. Such conclusions should be based on the 
results of analyses in a realistic design environment. 


8. COMPUTATIONAL ASPECTS OF THE FLUTTER TASK 


The purpose of this section is to delineate the computational aspects of 
the complete flutter task, which includes flutter analysis as well as struc- 
tural synthesis, i.e. the design of a structiore that satisfies the flutter 
requirements. Although the subject of this report is flutter optimization, 
i.e. structural synthesis aimed at a minimum weight structiire that satisfies 
the flutter requirements, it is useful to include flutter analysis, or flut- 
ter survey, in this discussion, since there is a large common data base and 
many common analytical tools . 

Structirral design aimed at satisfying flutter requirements must, to result 
in a viable airplane, also take into account strength requirements and require- 
ments related to manufacturing cost. Examination of a merit function that 
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combines structviral weight and manufacturing cost falls outside the scope of 
this study. The interface between structural synthesis with strength con- 
straints and synthesis with flutter constraints is discussed in Section T-5- 
There it is concluded that additional investigations in the areas of struc- 
tural modeling and structural optimization with combined flutter and strength 
constraints are needed before the best approach to total structiiral synthesis 
can be formulated. The organization of computer programs and modules dis- 
cussed in this section as a possible basis for definition of specifications 
for computer software envisions combining flutter optimization with satis- 
faction of stress constraints. 


8.1 The Complete Flutter Task 


The complete flutter task can be considered as being composed of three 
subtasks, subsequently to be discussed: 

1. Flutter survey of the original design. 

2. Initial structural resizing to satisfy all flutter requirements. 

3. Flutter optimization: weight minimization while explicitly satisfying 

flutter requirements and not violating strength requirements . 

8.1.1 Flutter Survey - The original design is defined by the external 
geometry, by a structural mass distribution derived by satisfying strength 
requirements, by an additional mass distribution representing fixed, non- 
structural airframe masses (e.g. powerplants, control system, furnishings) 
and mass distributions representing useful loads (e.g. fuel, payload). 

The flutter survey is a series of flutter analyses sufficient to identify 
any flutter deficiency that may exist in the original design over a range of 
operating conditions . Although details of the actual execution of the flutter 
survey may be different for different engineering facilities, it is believed 
that the flow diagram in Figure 8-1 is generally applicable. The general 
procedure is not new; it is discussed here in order to relate it to the overall 
design process. In Figure 8-1 oval boxes define engineer action points, 
although not necessarily manual operations; rectangular boxes represent com- 
puting modules. The computing process can proceed from one module to the next 
without engineer action, although an engineer's review may be inserted at 
any point . 

Starting with the design definition box in Figure 8-1, the engineer pro- 
ceeds to prepare, not necessarily manually, structural model data, inertia 
data and aerodynamics data. The structures module forms the stiffness and 
inertia matrices that are used in the vibration analysis. It performs, if 
necessary, the static coordinate reduction to keep the number of degrees of 
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Figure 8-1: Flutter Survey Task 








freedom for the vibration analysis within a practical limit. The structures 
module may also form the inertia matrix associated with the masses of the 
structural elements. If this is the case, the inertia data prepared by the 
engineer refer to nonstructural and useful load masses only. Inertia matrices 
for a number of useful load configurations, chosen on the basis of experience, 
and the output of the structures module are input into the vibration analysis 
module and vibration analyses, leading to natural frequencies and vibration 
modes, are performed. 

The natural frequencies and vibration modes can be reviewed by the 
engineer for checking purposes and, after the flutter analysis, for obtain- 
ing a better understanding of the flutter behavior of the airplane. In the 
case of a final design, the results of the vibration analysis can be com- 
pared with the results from ground vibration tests. From an analytical point 
of view, however, the vibration analysis is only necessary if the original 
niamber of degrees-of- freedom exceeds a practical limit for the flutter 
analysis. In that case, the vibration modes associated with the lower vibra- 
tion frequencies are, in general, used as generalized coordinates for the 
flutter analysis. The output of the vibration analysis can be formulated to 
include generalized stiffness and inertia matrices. In that case, the func- 
tion of the generalized-matrices module is to form only the generalized aero- 
dynamics matrices. However, to be more generally useful, this module should 
also include the capability of generating generalized mass and stiffness 
matrices by pre- and postmultiplication by modal matrices . This capability 
may be used if the generalized coordinates are not updated after each resizing 
step. In view of the options available for forming the generalized aero- 
dynamics matrices (Section 5.3), the generalized-matrices module may do more 
than a pre- and postmultiplication by the modal matrices obtained from the 
vibration analysis. Consequently, the output of the aerodynamics module 
may be a set of [hAw] matrices (Equation (5-l6)) for discrete values of the 
reduced frequency k, or may be a set of basic aerodynamic influence coeffi- 
cient matrices ^AIC(k)] and matrices [h], [dx] and [dz] ( Equations (5.8) 
and ( 5 . 32 )), or any combination between these extremes as discussed in 
Section ,h .k . The aerodynamics data preparation consists of defining aero- 
dynamics grid systems for the downwash collocation points and the aerodynamics 
loads points, and the selection of Mach numbers for which the flutter analysis 
will be performed. Structural grid data are input into the aerodynamics 
module to form the grid transformation matrices [h], [dx] and [dz]. 

static reduction of the stiffness matrix is a generally accepted method 
of reducing the number of degrees of freedom in the vibration analysis . 
Reference 36 presents an alternative that is worth, consideration. In the 
approach of Reference 36 , no reduction of the stiffness matrix takes place; 
in fact, the complete stiffness matrix is not assembled. Instead, the output 
of the structures module is a collection of submatrices that are used directly 
in the vibration analysis module to compute generalized stiffness and inertia 
data, and vibration modes. 
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From the generalized-matrices module the flutter survey procedure 
enters the flutter analysis module. The flutter analysis module should 
contain an interpolation routine for computing the generalized aerodynamics 
matrix for arbitrary k values . The flutter equation is solved by any 
suitable method and the output is a series of f-g-V diagrams (Figiire 3-1) 
for different Mach numbers, different distributions of useful load and for 
symmetric, anti-symmetric and possibly asymmetric modes. 

Figure 8-1 indicates three potential reanalysis loops . Any realistic 
procedure must account for the possibility that the initial choice of input 
parameters does not provide sufficient definition of the flutter character- 
istics of the original design. It is possible that the initial choice of 
k-values or speeds for which the flutter analysis is performed is insuffi- 
cient to determine flutter speeds or minimum damping in hump modes. Thus, 
review of the flutter analysis results may require return to the flutter 
analysis module for improved definition of the f-g-V diagrams. Review of 
the results may also Indicate the need for including more Mach numbers, 
inertia configurations associated with different useful load distributions, 
or stiffness matrices in the flutter survey. The possibility of Including 
more than one stiffness matrix follows from the fact that failed conditions 
must be considered. 

After sufficient flutter data have been generated, it is determined 
whether any flutter deficiencies exist. If there are no deficiencies, the 
flutter task is completed, unless design changes occur which make it neces- 
sary to repeat the flutter siirvey. Since the first survey resulted in 
engineering familarity with the flutter characteristics of the design, 
additional surveys usually can be restricted to fewer combinations of inertia 
configurations, Mach mmibers and failed conditions than were investigated 
during the first survey. 

If there are flutter deficiencies, a structural resizing is initiated 
which is aimed at satisfying all flutter requirements. 

8.1.2 Initial Structural Resizing - If the flutter survey of the original 
design indicates the presence of flutter deficiencies, it is desired to remove 
these deficiencies with a minimum weight penalty. Reference 1 indicates that 
it is theoretically possible to attain a minimum weight design by Judiciously 
adding structural mass in small quantities to those structural elements that, 
at each step, are most efficient in removing the deficiencies. For a struc- 
ture with a large number of design variables and several flutter deficiencies, 
this is, however, an impracticable approach. It is more efficient to first 
generate, in one or very few resizing steps, a structure without flutter 
deficiencies, but not necessarily with minimum weight, and then minimize the 
weight while avoiding flutter deficiencies. Most methods of flutter optimiza- 
tion discussed in this report are based on this approach. 

Two types of flutter deficiencies are recognized: l) too low a flutter 

speed (V^ < and 2) insufficient damping at the top of a hump mode 
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II 


(y, . y ,). Although both deficiencies are undesirable, the 

hump top > max allowed 

flutter speed deficiency occurs more frequently and is emphasized throughout 
this report. Optimization techniques aimed at satisfying flutter speed 
requirements can be generalized to include damping requirements . In this 
discussion, reference will be made, for convenience, to flutter speed require- 
ments, or constraints, only. 

Experience at the Lockheed-California Company, related to a realistic 
design enArironment , suggests that on the basis of engineering judgment, one 
flutter deficiency often can be identified as being most critical. That is, 
for a particular Mach number and useful load configuration there exists a 
deficiency , ' the removal of which is expected to result in all flutter defi- 
ciencies being removed. If this is not the case, then one or more additional 
applications of the following approach will lead to a design without flutter 
deficiencies. Neither case, in general, leads to an optimum design. 

Reference (l) indicates that a resizing column 
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( 8 . 1 ) 


where V is the most critical flutter speed, and C defines a magnitude 
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such that V = V_, is an efficient initial resizing. The column 
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is recognized as the gradient of the flutter speed. Other distributions 
of » however, may be considered, e.g.. 
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If it is difficult to define one most critical flutter mode, a resizing 
coliimn based on a weighted sum of two or more flutter speed gradients may be 
a better approach. In any case, initial resizing should take into account 
the efficiency with which design variables can increase the flutter speed 
and, thus, it is necessary, at this stage, to define design variables and to 
determine partial derivatives of the flutter speeds with respect to the design 
variables. In general, the initial resizing column can be defined as: 



(8.3) 


where the subscript j refers to the flutter speeds that are less than the 
required speed. 
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The initial resizing procediire, again, may be different for different 
engineering facilities and it probably depends, in its details, on the pre- 
ceding flutter s\xrvey procedure as well as on the subsequent optimization 
procedures to be used. With this in mind the essential features of the pro- 
cedure are shown in the flow diagram of Figure 8-2. Note that in Figure 8-2 
oval boxes still indicate engineer action points, but the rectangular boxes 
no longer define computational modules but computing activity in general, 
and no attempt is made to define specific modules. 


The endpoint of Figure 8-1 is the starting point for Figure 8-2: review 

of the flutter survey results. If there are flutter deficiencies, design 
variables may be defined and most critical flutter conditions selected. To 

av. 

generate the flutter speed derivatives characteristic roots and 

i 

vectors corresponding to the flutter points must be determined (e.g., by the 
two-dimensional Regula Falsi followed by a subroutine for detennining char- 
acteristic vectors). In addition, derivatives of the mass, stiffness and 
aerodynamics matrices are required. If there is a static reduction of the 
stiffness matrix, the static reduction matrix must be explicitly generated 
in order to compute the derivatives of the reduced stiffness matrix. The 

av. 

flutter speed derivatives can be computed following the formulation of 

i 

Reference l4, a compact version of which is included in Reference 1. The 

av. 

analyst may want to review the values of before deciding on the direc- 

i 


tion of the resizing column, 


ri 

lam. 
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, or it may be formed automatically 


by the program. Incremented Flutter Analysis (References 1 and 4) , is a con- 


venient method of determining each C. that results in V. = V . If 

j j Rj 


/'av.' 

j 

\dm. 
mines 


contains only positive elements, the largest value of C. deter- 
1 /J r 1 

the resizing column according to equation (8.2) which results 


in all V. > V . If there is uncertainty about this result, the structure 
J R-i 


is incremented by 


Ami) 


and a new vibration and flutter analysis is per- 


formed for selected combinations of Mach number and useful load. If necessary. 


more critical flutter speeds V^ are selected and the process is repeated. 

The procedures described, and illustrated in Figure 8-2, are aimed at a 
one-step resizing to reach the goal of satisfying all flutter requirements. 
Variations of this procedure are possible but would include the same basic 
computational modules . 

•• 
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8.1.3 Flutter Optimization - Of the three suhtasks discussed in this section, 
the flutter optimization task is most dependent on the computational methods 
adopted by an engineering facility. However, this method dependency is 
mainly concentrated in the actual resizing procedure with a constant set of 
generalized coordinates (invariant vibration modes). The structural analysis 
and vibration analysis, to be repeated several times during the optimization 
task, can be defined in general terms, relatively independent of the resizing 
procedure used. 

The methods of optimization discussed in Section 6 deal with only a 
small aspect of the optimization task, namely, the resizing when given a con- 
stant set of generalized coordinates i ( i. e. , invariant vibration modes). It is 
most likely that modal updating is necessary (Section k) . If a static reduc- 
tion of the stiffness matrix is used to reduce the number of degrees of 
freedom in the vibration analysis (Section T-l), the reduced stiffness matrix 
may be a nonlinear function of the design variables . This results in addi- 
tional computations to determine derivatives of the stiffness matrix and 
possibly the inertia matrix. 

In this section the overall flutter optimization task, except for the 
resizing for a constant set of generalized coordinates, is delineated under 
the following assumptions: 

1. There will be modal updating. 

2. Static reduction of the stiffness matrix is required. 

3. For each resizing step the derivative of the reduced stiffness matrix is 
determined exactly (Equation (7-6)). 

U. No static reduction of the mass matrix is required. 

5. There is one active flutter speed constraint. 

6. Strength requirements are satisfied. 

Keeping in mind these assumptions, a flow diagram is formulated (Fig- 
ure 8-3) that delineates those computational steps that are considered 
independent of the method chosen for determining the resizing steps, given 
the flutter speed derivatives. Several of the computations indicated in 
Figure 8-3 are identical to those in Figure 8-2. 

The inputs into the flutter optimization task are: the starting stiff- 

ness and mass matrices that are output from the initial resizing task; the 
incremental stiffness and mass matrices; the basic aerodynamics input; and 
the aerodynamics derivatives input. 

The starting (first current) stiffness and mass matrices are used to 
define generalized coordinates via a vibration analysis. These genera]. ized 
coordinates along with the associated generalized stiffness and mass matrices. 
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Figure 8-3: The Flutter Optimization Task 

















and the basic aerodynamics input are used to obtain a point solution of the 
flutter equation and the associated characteristic vectors . The character- 
istic vectors are combined in turn with the incremental stiffness matrix 
pre- and postmultiplied by the static reduction matrix [gr], with the 
incremental mass matrix, and with the aerodynamics derivatives input, to form 
the derivative scalars (See Reference l). The flutter speed derivatives are 
formed and input into the resizing module. Other inputs into the resizing 
module depend on the method of optimization used, but include minimum size 
constraints equal to the original design, i.e. , before the initial resizing, 
and program control parameters to be chosen by the analyst. 

The resizing module generates a column of design variable increments 
defining one step in a resizing process that usually comprises- several steps. 
The total mass associated with the new values of the design variables is 
compared with the previous total mass. If the total mass has not converged 
to a minimum, the design variable increments are used to generate new current 
stiffness and mass matrices and the process is repeated. 

Satisfying the stress constraints can be accomplished in various ways. 

One approach, which does not involve the resizing module, is shown in 
Figure 8-3. In it, when a minimum total mass is reached, the loads and 
stress analysis is redone. Due to stiffening of the original design, there 
may be stress violations due to redistribution of internal loads or changes 
in the external loads. Where such violations occur, the element sizes are 
increased to satisfy the stress constraints and new current stiffness and 
mass matrices are formed and the flutter optimization is repeated using 
updated minimum sizes. After a new convergence on a minimum total mass, the 
stresses are checked again and, if necessary, the entire process is repeated. 
The increased element size will, in general, cause a change in flutter speed. 
This is expected to be small, since the stress violations are most likely to 
be in elements that have not increased in size during the flutter optimiza- 
tion and, thus, are ineffective in changing the flutter speed. 

In another approach there is a strong interaction between a stress 
analysis module and the resizing module such that each resizing step is con- 
strained such that stress constraints are not violated. Two possibilities 
can be distinguished: the loads are assumed constant or the loads are recal- 

culated at each step. In the former case the loads and stress analysis shown 
in Figure 8-3 must follow the flutter optimization as described in the previous 
paragraph . 

It is noted that the aerodynamics input, in Figure 8-3, is not defined 
in terms of specific matrices. The discussions in Section 5 and Reference 1 
indicate that there are many options for formulating the matrices of general- 
ized aerodynamic force coefficients, and it is considered outside the scope 
of this discussion to present a definite choice. It is worthwhile noting, 
however, that if cubic spline interpolation is used for the aerodynamics, 
the basic aerodynamics input and the aerodynamics derivatives input are 
identical matrices. 



Since no recommendations for a specific method of optimization are made 
as a result of this study, a further definition of the resizing module, in 
Figure 8-3, is considered to fall outside the scope of this report. 


8.2 Aspects of the Computing System 


General aspects of the computing system required to perform the flutter 
task are discussed without attempting to define detailed specifications for 
such a system. 

When properly divided into well defined building blocks the complete 
flutter task, including the optimization, is relatively straightforward, 
inviting extensive automation. The flutter optimization task seems especially 
well suited for complete automation. It is questionable, however, whether 
the actual design experience available is sufficient to decide on all aspects 
of the optimization task and to embark on the design of a computing system 
that can handle efficiently a structural design that is defined by several 
thousands of finite elements . In view of this the Lockheed-California Company 
has first developed a semi-automatic system, based on its Computer Graphics 
system. It has been used on the configuration upon which the numerical 
examples of Appendix A are based, on an arrow wing supersonic transport study 
(Reference l) and on an actual hardware problem. The authors believe, how- 
ever, that a batch process system with maximum automation options is a 
desirable design asset, even if in its initial version it is somewhat 
restricted in the number of structural degrees of freedom, design variables 
and flutter constraints it can handle. In this section some aspects of a 
batch process computing system are discussed. 

An engineering facility that already has a computing system for flutter 
analysis (flutter survey task) may want to restrict its batch process system 
to the initial resizing task and the flutter optimization task. This, how- 
ever, seems only justified if the data input and output of the existing flut- 
ter analysis system can easily be made compatible with the input requirements 
for the other tasks. Since the repetitive analyses associated with flutter 
optimization put extra emphasis on computational efficiency, a facility may 
decide to update its flutter analysis system as part of the introduction of 
a flutter optimization capability. In the following it is assumed that a 
batch process system for the complete flutter task is to be designed. 

When comparing the flow diagrams in the Figures 8-1 through 8-3 it is 
clear that the three tasks represented in these three figures have a large 
number of computational functions in common. Thus, the first point of con- 
sideration is whether three independent programs should be developed within 
an existing computing system or whether one program, or a new system should 
be developed for the complete flutter task. The choice depends on what 
computing system is available at a facility and, to a certain extent, personal 
preferences of the engineers. From a practical engineering point of view it 
seems self-evident that, whatever course is chosen, data format compatibility 
is mandatory and that a data management system is required. 
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In selecting a computing system for the complete flutter task, con- 
sideration should he given to a system that has access to an existing matrix 
algebra computing system. If this access is not available, the complete flut- 
ter task system should include some generalized matrix algebra capability to 
enable the user to depart from a rigid format. 

In designing a computer system for the complete flutter task, two 
approaches can be distinguished: Self-Contained Programs and Variable Job 

Stepping. 

In the first approach, which is the more common of the two, there may 
be a self-contained program for each of the subtasks comprising the complete 
flutter task, or some or all subtasks may be combined into one self-contained 
program. Each program has its own executive module which controls calling 
into core the various computing modules (sub-programs) as they are needed 
dturing the entire computer operation for that program. Organizing the com- 
plete flutter task in one program with one executive module would result in 
a very large program with some complex input/output interface problems as well 
as core overlay problems . If existing batch programs are to be integrated 
into such an overall computing program, extensive modification ml^t be needed 
in order to resolve some of these problems. If the complete flutter task is 
covered by more than one program, automatic transfer from one program to 
another might prove impracticable, thus limiting the overall automation 
attainable. Since engineer reviews are required during the computational 
effort, however, this may not be an important limitation. Consideration must 
be given to the degree of commonality between modules in the separate pro- 
grams that perform the same function. Failing to achieve complete commonality, 
e.g., due to different overlay requirements, increases the effort needed to 
update a functional module. 

Variable Job Stepping consists of a number of separate computer programs 
each representing a computing module, such as those defined in Figure 8-1, 
which are controlled by another program called the Executive. Variable Job 
Stepping, therefore, is a sequence of separate computer job steps in which 
each job step is a function of the preceding job step(s) as determined by 
the Executive Module. The Executive Module (a separate program) monitors 
the task completion codes of the preceding steps and, based on the instruction 
code supplied by the engineer, determines which computing module is next 
required. Prior to transferring the control to the computing module, the 
Executive Module prepares the input data for that computing module in accor- 
dance with its data format requirements. The Variable Job Stepping system 
is inherently modular in approach. To add another module, only the Executive 
Module (program) needs to be modified and reloaded as a new executable pro- 
gram in addition to the new computing module. Existing batch process programs 
can be included in a Variable Job Stepping system with little modification to 
the batch programs. If a new method for computing, say, aerodynamic matrices, 
becomes available, again only the Executive Module needs to be modified. It 
is also worth noting that each program within the Variable Job Stepping system 
could be executed as an independent program outside the Executive Module 
Control. Figure 8-U illustrates the principle. The Executive is loaded and 
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Figure 8-4; Computing System Using Variable Job Stepping 





calls in the next module needed, say C, which is loaded over the same space 
as was occupied by the Executive Module. When module C has completed its 
task it calls back the Executive Module, at the same time giving it instruc- 
tions that depend on the output of module C. The Executive Module calls in 
the next computational module, again while annihilating itself from core. 

Which of the two approaches discussed is preferred depends on many 
factors and, without a more in-depth look at the details of the computing 
system, cannot be determined at this time. 

As stated above, the complete flutter task is relatively straightfor- 
ward, when properly defined, and will, in principle at least, not give rise' 
to great programming difficulties. The great challenge in designing the 
system lies in maximizing its capacity, in terms of number of structiiral 
coordinates, design variables and flutter constraints, while keeping computing 
cost reasonable. This does require a very efficient use of all computer 
resources: maximum utilization of available core, efficient input/output 

routines , systematic labelling and external storage of data blocks , efficient 
compacting of data blocks (e.g. sparse matrices), etc. Obviously, close 
cooperation between exp.erts in flutter analysis, optimization procedures and 
computer programming is required in order to obtain an efficient computing 
system. 


9 . CONCLUSIONS AND RECOMMENDATIONS 

Based on the investigations and evaluations performed during this study 
and augmented by the additional supporting activities conducted concurrently 
at the Lockheed-California Company, it is concluded that an efficient, 
practicable flutter optimization module can be formulated and implemented 
with a reasonable amount of further development. Some of the particular con- 
clusions supporting this general conclusion are presented in the following: 

Aerodynamics parameters may be efficiently represented by the five- 
matrix product shown in equation (5.12). The most efficient method of imple- 
mentation of this equation depends upon a complex relationship involving the 
numbers of aerodynamic integration points, downwash points, generalized modal 
coordinates, discrete structural coordinates and reduced frequencies used, 
as well as the interpolation procedures employed. To be completely general, 
several options should be available in order to provide alternate procedures 
for the generation of these parameters . In relation to the complete flutter 
optimization task, however, it is concluded that severely limiting the number 
of such options would not significantly increase the required computing 
resources . 

Repetitive flutter solutions of the type required in structural resizing 
procedures may be efficiently obtained using the 2-D Regula Falsi method. 
Although other flutter solution procedures might be developed to the same 
degree of efficiency and reliability demonstrated by the Regula Falsi in 
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obtaining repetitive flutter solutions, it is concluded that this procedure 
has a wider range of application than other procedures considered. The 2-D 
Regula Falsi method can be used in the solution of a general fomi of the 
flutter equation in order to obtain the values of any two dependent variables 
required to satisfy the equation. It can be used in a direct form of Incre- 
mented Flutter Analysis to solve for the magnitude of the increment of a 
specified design variable, along with an additional variable, necessary to 
satisfy a flutter constraint. Other applications, such as the determination 
of the minimum damping of a hump mode, are also possible. 

Resizing procedures vary widely in particular detail, with each proce- 
dure considered exhibiting one or more advantageous features. In terms of 
a realistic ; design effort , however, it is concluded that the general differ- 
ences between the arbitrary step-size procedures and the defined step-size 
procedures are more significant than the differences between the individual 
procedures in each category. Based on the numerical evaluations of the 
idealized test case of Appendix A, it would appear that the arbitrary step- 
size procedures produce the same mass reduction as the defined step size 
procedures with less total computing cost. It is not clear that this same 
result would be obtained for a more complex design effort. It is concluded, 
however, that either type of procedure can be included in a flutter optimiza- 
tion module with no undue difficulty. 

Practical considerations in the implementation of a flutter optimization 
procedure can have at least as profound an effect on the performance of the 
flutter module as does the selection of the three major elements considered 
thus far. These considerations include the choice of structural model, num- 
ber of structural degrees of freedom retained, method of generation of 
incremental stiffness matrices, number of modal coordinates used, frequency 
of updating vibration modes and frequency of updating flutter derivatives . 

No definite conclusions are available regarding these considerations, since 
it was not possible to conduct the required Investigations within the scope 
of the present study. 

The conclusions reached in the course of the present study, and presented 
above, lead to the following recommendations for development of a flutter 
optimi zation module : 

1. Computer coding for the aerodynamics submodule should be based on the 
five-matrix product of equation (5.12). The number of options in form- 
ing the five-matrix product should be limited to those forms that are 
most generally useful. 

2. The flutter solution submodule for the repetitive flutter solution pro- 
cedure should be based on the two-dimensional Regula Falsi approach. 

A global solution procedure, such as the p-k method of Reference 3 or 
the Desmarais -Bennett method (Reference 7) should be included for the 
initial and final flutter survey. 
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3. 


At least two candidate resizing procedures should be evaluated in a 
realistic design effort such as the arrow wing flutter optimization task 
performed under Contract NAS-1-12288. One resizing procedure should be 
of the arbitrary step-size type and the other a defined step-size type. 

To facilitate the comparison, it would be useful to maintain as much 
commonality between the two methods as possible. As an example, the 
resizing procedure formulated in Section 6.J.h could be used as the 
defined step-size procedure, and the gradient projection search of 
Rudisill-Bhatia (Section 6.2), modified to incorporate the move vector 
of Section A. 3. 3, could be used as the arbitrary step-size procedure. 

U. Depending on the outcome of these evaluations, specifications should 

be developed for the selected procedure and the required computer coding 
accomplished. It should be recognized that it may be desirable to retain 
the option of using either type of resizing procedure in the final 
module . 

5. Investigation of the flutter optimization task should be extended to 
include those problems encountered in a realistic design environment 
which have a direct influence on the performance of an optimization 
procedure. These problems are discussed in Section 7 and mentioned 
briefly in the conclusions above. To be of most use, it is felt that 
such investigations must be made using a structural design task of the 
complexity of the arrow wing study performed by Lockheed for NASA 
(NAS-1-12288) . 

6. The extent to which it is feasible and desirable to integrate the flutter 
optimization task with the strength optimization task should be investi- 
gated with emphasis on means of simplifying the formulation of the 
strength constraints . 


APPENDIX A 

NUMERICAL EXAMPLES OF RESIZING PROCEDURES 


A.l INTRODUCTION 


To provide the basis for a direct comparison of the several candidate 
resizing procedures (Section 6) in performing a simplified flutter optimiza- 
tion task, an idealized test case was formulated and numerical evaluations 
were conducted. Not all candidate procedures were evaluated to the same 
degree. In most cases the process was discontinued as soon as the results of | 
interest were obtained. 
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The structural model on which the test case was based is a simple El, 

GJ beam representation of a subsonic transport airplane. Although no attempt 
was made to simulate a realistic design process in detail, the resulting 
optimization might be regarded as typical of the preliminary design phase of 
the development of such an airplane. The test case was designed to avoid the 
difficulties associated with modalization (Section 4) and the nonlinear 
stiffness effects (Section T.l) encountered in more practical optimization 
efforts . 

In implementing the various resizing procedures, no specialized computer 
programs were formed. Instead, existing batch, graphics and remote terminal 
systems and programs were employed, augmented by hand computations where 
necessary. As a result, no direct comparison of the computer resources 
,req.uired by the various resizing procedures is available; such information 
relative to this idealized test case is of little practical significance in 
any event. 


A. 2 STRUCTURAL MODEL 


The structural model used for the test case is an El, GJ beam representa- 
tion of a subsonic transport airplane. 

There are 9 grid-points on the wing semi-span, 12 grid-points along the 
fuselage center line (Figure A-l) and other miscellaneous grid-points used to 
define a rigid empennage and control surface, for a total of 6j elastic 
degrees of freedom. Symmetric boundary conditions are imposed. All inertia 
and aerodynamics coordinates are retained. Eighteen elastic degrees of 
freedom, not associated with design variables, are eliminated by static stiff- 
ness condensation, reducing the number of degrees of freedom to 49 . 


The design variables are the torsional stiffnesses of the eight struc- 
tural elements indicated in Figure A-l. The bending stiffness of the wing is 
not varied. The design variables are defined as increments over values 
defining a simulated strength design and the associated inertia and stiffness 
matrices are expressed in terms of a unit mass , so that the total inertia and 
stiffness matrices are as expressed in equations (A.l) and (A. 2). 



(A.l) 


(A. 2) 
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Figure A-1: Structural Model 

are the matrices of the fixed stiffness and inertia 
are the stiffness and inertia matrices associated 
with a unit mass of the design variable m^ . For the present piirposes, the 
minimum-size or strength designed portions of the design variables are 
included in the K^, matrices, so that the AK^ , AM^ matrices include 

only the design variable increments relative to the strength design. 

The strength designed configuration is a fictitious configuration defined 
such that a flutter problem is assured within the design envelope of the air- 
plane. The characteristics of this critical flutter root are shown in Fig- 
ure A-2, indicating a flutter speed of approximately 226. m/s EAS (khO KEAS) . 
Using this as a starting point, the test case required that the flutter speed 
be increased to 270.1 m/s EAS (525 KEAS) with a minimum of weight increase in 
the torsional stiffness design variables. 


where 

and 


H 


and 

and 


w 

M 
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Figvire A-2*. Frequency and Dan 5 )ing of Critical Flutter Root 


In addition to the flutter-deficient (226.4 m/s EAS) configuration, two 
auxiliary configurations are required. For the constant flutter speed pro- 
cediores, a configuration having the required flutter speed (2T0.1 m/s EAS) 
hut a non-optimijm distribution of the design variables is needed. This was 
obtained by increasing the design variables in a manner equivalent to raising 
the torsional stiffness of the wing by a uniform factor until the required 
flutter speed was reached. For the penalty function procedure (Reference l6) , 
an initial configuration having a flutter speed in excess of the final flutter 
speed is required. This configuration was obtained in the same manner as the 
previous one, except that the flutter speed was increased to 280.4 m/s EAS 
(545 KEAS). The design variable distributions and total design variable mass 
for these configurations are shown in Table A-1. A man- in- th e-loop resizing 
procedure, different from any of the methods discussed in Section 6 and using 
Incremented Flutter Analysis as the principal tool, indicated that the mini- 
miim mass for the required flutter speed is approximately 24p.6 kg (550.2 lbs). 
This value is used as a bench mark for further comparisons . 
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Flutter 



Design Variable Mass, kg* 



Speed 
m/s EAS 

1 

2 

3 

4 

5 

6 

7 

8 

Total 

226.4 

0 

0 

0 

0 

0 

0 

0 

0 

0 

2T0.1 

113.5 

105.2 

90.4 

65.5 

63.6 

31.3 

21.8 

15.3 

506.6 

280.4 

141.3 

131.0 

112.5 

81.6 

79.2 

38.9 

27.1 

19.1 

630.9 ■ 


* In accordance with the definition of design variable this is a mass 
increase above the simulated strength level design. 


Table A-1: Initial Configviration; Non-Optimum Distribution 


A. 3 METHOD OF SIMDDYNES 


The method of Simodynes is reported in Reference 15 and discussed in 
Section 6.3 of this report. The numerical evaluations reported here are in 
fact evaluations, or partial evaluations, of three distinct methods. The 
Simodynes method itself was evaluated only to the extent of the first 
resizing cycle. At that point, two undesirable features of the method were 
identified and a modification of the method was implemented. The numerical 
evaluation was then continued, using this modified method. In the modified 
method one of the objectionable features of the original method, the frequency 
constraint, is avoided, as is discussed in Section 6.3. Finally, a second 
modification was incorporated by which the influence of the choice of the 
dependent design variable on the resizing step is eliminated, and an addi- 
tional numerical evaluation was performed lising this procedure. 

A. 3.1 Original Method of Simodynes - As discussed in Section 6.3j the method 
of Simodynes is a procedure for minimizing the total weight of a set of 
design variables while maintaining a fixed flutter speed and frequency. In 
order to satisfy these constraints, two of the design variable masses are 
considered to be dependent functions of the remaining design variable masses. 
The resizing direction is then determined by the gradient of the total weight 
subject to the flutter speed and frequency constraints . 

An initial resizing step was generated for an arbitrarily chosen total 
mass reduction, of kg (lOO lbs) for each of three pairs of depen- 
dent design variables: 1 and 2, U and 5» and T and 8. The resulting values 

of the design variable increments are shown in Table A-2. 































I 


Dependent 
Design Variables 

1 

D 

2 

“sign 

3 

/■ariable 

4 

Incre 

5 

ments , 
6 

kg 

7 

8 

Total 

1 and 2 

-198.6 

155.9 

-0.0 

0.0 

-0.0 

- 0.3 

-1.0 

- 1.4 

-45.4 

U and 5 

-0,2 

-0.2 

-0.2 

-106.4 

68.3 

- 0.5 

- 2.3 

- 3.9 

-45.4 

7 and 8 

-10.1 

-10.0 

-9.7 

-9.0 

-7.6 

- 4.3 

20.6 

-15.3 

- 45.4 


Table A- 2 : Initial Design Variable Increments (first resizing step) 

for Three Pairs of Dependent Design Variables 


The dual constraint of flutter speed and frequency leads to relatively 
large positive and negative increments in the dependent design variables. 

When design variable pairs 1 , 2 and 4 , 5 si’s chosen as dependent design 
variables the magnitude of the negative increment is larger than the avail- 
able amount (cf. Table A-l) . The question arises whether to invoke the 
minimum size constraint and accept the resulting larger drift in flutter 
speed, or recognize the violation of the minimum size constraint as a tem- 
porary situation that will be corrected in subsequent steps, or reduce W^. 

It should be noted, however, that as a result of the large increments of the 
dependent design variables, the amo^^nt of resizing towards the goal of minimum 
total mass occurring in the independent variables is small when the pairs 1, 2 
and U, 5 are the dependent design variables. How serious a drawback is 
implied by these considerations has not been pursued in detail, since it is 
believed that the frequency constraint by itself puts Simodynes method at a 
distinct disadvantage, certainly in view of the rather straightforward manner 
in which this constraint can be removed, as is demonstrated in the modified 
Simodynes method which is the subject of the next section. 

A. 3. 2 Modified Simodynes Method - The method of Simodynes was modified to 
eliminate the flutter frequency constraint, substituting flutter frequency 
as a dependent variable in place of one of the two dependent design variables 
used in the original method (Section 6 . 3 ). , 

As in the case of the original method, an initial resizing step was 

generated for a total mass reduction of = ^ 5.4 kg (lOO lbs) for each 

of three choices of the dependent design variable. The results, shown in 
Table A- 3 , indicate more reasonable distributions of the resizing increments 
due to the elimination of the frequency constraint, although the sensitivity 
to the choice of dependent design variable remains. Specifically, the restilt 

9 V 

of choosing as a dependent design variable one for which ^ — is small in 

dm 
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Dependent 



Design 

Variable Increments 

, kg 



Design Variable 

1 

2 

3 

k 

5 

6 

7 

8 

Total 

1 

-50.1 

0.0 

0.1 

0.4 

0.9 

1.4 

1.4 

0.4 

-45.4 

5 

-12.8 

-12.0 

-10.8 

-6.1 

-14.2 

8.5 

8.8 

-6.8 

-45.4 

8 

-3.k 

-3.0 

-2.3 

0.4 

3.9 

8.8 

9.0 

- 58.8 

-45.4 


Table A-3t Initial Design Variable Increment (first resizing step) 
for Three Dependent Design Variables; Modified 
Simodynes Method 


magnitude is demonstrated in Table A-3- With number 1 as the dependent 
variable, a large negative increment of that variable is required to balance 
rather small increments in the independent design variables in keeping the 

flutter speed constant on a linear basis. For a given weight decrement W^, 

the amount of resizing towards the goal of minimum total mass occurring in 
the independent variables is small. It would seem, therefore, that it is 
important to choose as a dependent design variable one that corresponds to 

0V 

a large value With the proper choice of this variable, it would appear 

that a reasonably efficient procedure might result. In order to assess this, 
additional resizing steps were executed, using design variable 5 as the 
independent variable. The flutter speeds and design variable distributions 
for six resizing steps are presented in Table A-k . One other minor departure 
from the original method is that the flutter speed was adjusted to approximately 
270. 1 m/s EAS at the beginning of each resizing cycle by a uniform percentage 
increase in the current design variable distribution. 

It is believed that the efficient performance indicated by the results 
shown in Table A-k may be somewhat better than might typically be achieved 
by this method. These evaluations were performed relatively late in the con- 
tract effort, so that much experience with this test case had been accumulated. 
This allowed a choice of total mass reduction steps which minimized the 
required number of steps to achieve the mass reductions shown. Specifically, 

the successive values of W^ chosen (90. Tj 90.7* 90.7, ^5*^, 9-1, 9.1 kg) 

anticipate the minimum weight. 

A. 3. 3 Improved Move Vector - To eliminate the sensitivity of the modified 
method to the choice of dependent design variable, a modified move vector 
was fonnulated and tested in a further evaluation. 
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Step 

1 

2 

3 

Des: 

U 

Lgn Va] 

5 

'■iable 

6 

Mass , 
7 

CO 

Total 

Flutter Speed 
m/s EAS 

0 

113.5 

105.2 

90. U 

65.5 

63.6 

31.3 

21.8 

15.3 

506.6 

270.1 

1 

88.0 

8l.l 

68.8 

53.1+ 

35.2 

1+8.3 

39.3 

1.7 

415.9 

267.2 

2 

6o.T 

55.1 

UU.T 

38.5 

71.1+ 

1+8.9 

35 . 1 + 

0.0 

354.7 

268.3 

3 

2k.6 

22.2 

17.2 

30.5 

57.2 

63.0 

1+8.1 

15.9 

278.9 

267.8 

k 

6.1 

5.5 

3.6 

26.3 

89.8 

6 I +.5 

47.2 

6.1 

249.1 

268.6 

5 

0.6 

1.0 

O.k 

27.5 

93.8 

68.1 

1 + 9.8 

7.6 

248.8 

269.7 

6 

0.0 

0.0 

0.0 

28.2 

91.3 

69.9 

51.2 

8.6 

249.2 

270.0 


Table A-h: Mass Reductions with Modified Simodynes Proceditre 
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The new move vector is shown in eq.uation (A.3)» where -r — is the 

dm. 

1 

derivative of the flutter speed with respect to the design variable m^ . 

The first, positive, component of this move vector will be recognized as the 

(A. 3) 


’ 

Am, 

1 


av 

. — "h . 

1 

■“ 3.' 

am. 

1 


av/am^ 




a = 


b = 


l^ij |av/am^| 

aj^av/am^ J_hH 




(A.U) 


(A. 5) 


velocity gradient. The scalar "a" is determined from equation (A.U) , where 
is 

{Ami} • 


is the total mass specified for the positive component of the vector 


The scalar b defines the magnitude of the negative component of 
the vector such that the change in flutter speed, on a linearized 
basis, is zero. 


10 T 


such that the change in flutter speed, on a linearized 



The flutter velocity derivatives, needed for this move vector, are not 
computed directly in either the Simodynes or modified Simodynes procedure. 
Instead, partial derivatives of the dependent variables with respect to 
independent variables are generated. For the modified Simodynes procedure, 

considering the dependent design varialbe, m^, and one independent design 

variable, m^ , keeping all other independent variables constant, the condition 

V(m , m. ) = constant leads to: 
u 1 

8m 8V/8m. 

u _ 1 _ 8v 

8m. 8V/8m 8m. 

1 u 1 

8V 

where is a normalized flutter speed derivative, 

i 

complete set of such derivatives is then equivalent to the flutter velocity 
derivatives normalized to the flutter velocity derivative of the dependent 
design variable. Examination of equations (A. 3), (A.k) and (A. 5) shows that 

the column of resizing increments, Am^ , remains invariant with normalization 

of the flutter velocity derivatives. Accordingly, the normalized derivatives 
indicated in equation (A. 6) are used to evaluate the resizing increments 
shown in equation (A.3). 

Using the move vector described here, with a value of = 45.^ kg 

(lOO lbs), the results of eight resizing steps are shown in Table A-5. In 
comparing these results with the results in Table A-4 and noting that the 
"improved" move vector requires more steps to reach the minimum total mass , 
reference is again made to the fact that the results of Table A-4 were gen- 
erated after considerable experience had been gained with this idealized test 

case. Specifically, it was found that the value of chosen affects the 

convergence of the procedure. Note that the quantity in the present 

procedure has no direct relationship to the quantity used in the modified 

Simodynes procedure . 

To determine the effect of varying W^, steps 5> 6 and 7 were repeated 
with values of of 90.7 j 90.7 and 136 kg, respectively. The results, 

presented in Table A-6, demonstrate that improved performance of the method 
co\ild be ejq^ected for a more favorable choice of the value of W^. 


(A.6) 

The negative of the 
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Design Variable Mass, kg 


Flutter Speed 

Step 

1 

2 

3 

4 

5 

6 

7 

8 

Total 

m/s EAS 

0 

113.5 

105.2 

90.4 

65.5 

63.6 

31.3 

21.8 

15.3 

506.6 

270.1 

1 

58.7 

67.5 

66 .3 

58.6 

64.8 

39.4 

30.1 

6.9 

392.3 

269.2 

2 

25.0 

43.2 

50.7 

55.4 

68 .0 

46.8 

36.6 

6.8 

332.8 

269.6 

■3 

0.3 

24.6 

38.5 

53.2 

71.4 

53.1 

4 i .7 

7.9 

290.7 

269.8 

1; 

0.0 

7.5 

26.9 

50.8 

74.0 

58.0 

45.3 

8.5 

271.1 

269.9 

5 

0.0 

0.0 

i 6.5 

48.6 

76.2 

61.9 

47.9 

8.8 

260.0 

270.0 

6 

0.0 

0.0 

7.0 

46.4 

78.1 

64.8 

49.6 

8.9 

254.8 

270.0 

T 

0.0 

0.0 

0.0 

44.5 

79.6 

67.2 

50.9 

9.1 

251.2 

270.0 

8 

0.0 

0.0 

0.0 

42.1 

80.4 

68.5 

51.3 

8.6 

250.8 

270.1 


Table A -5 •' Mass Reductions with Improved Move Vector 






Design 

Variable Mass, kg 



Flutter Speed 

Step 

kg 

1 

2 

3 

4 

5 

6 

7 

8 

Total 

m/ s EAS 

5 

90.7 

0.0 

0.0 

5.0 

45.6 

77.7 

65.1 

49.8 

8.7 

252.0 

269.8 

6 

90.7 

0,0 

0.0 

0 .0 

40.9 

80 .2 

69.0 

51.5 

8.3 

250.0 

270 .0 

7 

136.1 

0.0 

0.0 

0.0 

35.1 

82.2 

71.9 

52.0 

1 

8.4 

249.6 

270.0 


Table A-6: Mass Reductions with Increased Values of 


A.k GRADIENT METHODS OF RUDISILL AND BHATIA 

In Reference l 4 , Rudisill and Bhatia present a method of generating 
flutter velocity derivatives and they suggest several resizing procedures 
using these derivatives. These procedures are discussed in Section 6 . 2 . Two 
of these procedures, the velocity gradient search and the gradient projection 
search, were chosen for nTomerical evaluation and the results are presented in 
the following sections. 


109 



A. 4.1 Velocity Gradient Search - The velocity gradient search is used to 
increase flutter speed by a series of resizing steps in which the increments 
in the design variables are proportional to the corresponding elements of the 
velocity gradient. The resulting distribution is not optimiom, but is a good 
initial distribution for procedures in which the total mass is minimized at 
constant flutter speed. One obvious application of this procedure is in 
increasing the flutter speed of a flutter-deficient design to the required 
flutter speed. Starting with the 226.4 m/s EAS configuration of the test case 
structural model, the flutter speed was increased to approximately 270.1 m/s 
EAS in five steps (Table A- 7 ). The design variable increments were formed 
according to equation {A.j) , where the nominal velocity increments, AV, were 


Am. 

1 


AV 

av/am^ . av/am^ ' 


av 

am. 

1 


(A. 7 ) 


12 . 9 , 10 . 3 , 10 . 3 , 5.1 and 3.1 m/s EAS. The actual velocity increments do 
not correspond exactly to the nominal increments because of l) the nonlinear 
relationship between the increment in flutter speed and increments in the 
design variables, and 2) the fact that the velocity derivatives, as defined 
in Reference l 4 , are not based on matched atmospheric conditions of Mach 
number, speed and altitude, whereas the flutter equation was solved directly 
for matched conditions . 


To evaluate the effect of the number of resizing steps used to produce 
a given velocity increment, the magnitude of the original increment was 
determined which resulted in the same flutter speed of 270.5 m/s EAS in one 
resizing step. The required mass, 292.3 kg, does not differ greatly from 
the multi-step result of 286.7 kg, when compared with the optimum total 

design variable mass of , 249.6 kg pounds (for V^ = 270.1 m/s EAS). 






Design Variable ffess 

. kg 


Flutter Speed 

Step 

1 

2 

3 

4 

5 

6 

7 

8 

Total 

m/ s EAS 

1 

1.1 

1.6 

2.4 

5.7 

10.1 

16.6 

17.8 

8.1 

63.4 

239.0 

2 

2.6 

3.9 

6.0 

13.6 

23.2 

33.2 

32.6 

15.7 

130.9 

249.9 

3 

5.0 

7.4 

11.2 

24.6 

40.3 

51.3 

47.7 

23.6 

211.0 

261.1 

4 

6.7 

9.8 

i 4.7 

31.5 

50.4 

60.8 

55.4 

27.5 

250.8 

260.9 

5 

7.9 

11.6 

17.2 

36.2 

57.1 

66.8 

60.1 

29.8 

286.7 

270.5 


Table A- 7 ‘ Design Variable Mass and Flutter 
Speed; Velocity Gradient Search 
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k.h.2 Gradient Projection Search - The gradient projection search is a 
procedure for reducing the total mass of the design variables while attempting 
to maintain a constant flutter speed. The column of resizing increments is 
derived from the velocity gradient and the mass gradient as indicated by 
equation (A. 8) 



8M 


9m. 

1 


+ X. 


9V 


9m. 

1 


|9V/9mJ |9M/9m^| 
|9V/9m^J |9V /9m^J 


(A. 8) 


(A. 9) 


where the scalar 


is determined as in equation (A. 9). Comparison with 


equations (A. 3), (A.U) and (A. 5) shows that this resizing column is similar 
to that of Section A. 3. 3 except the mass gradient is used in place of the 
reciprocals of the flutter derivatives . For comparison with this previous 


procedure, the product coefficient of the velocity gradient, was 

chosen to give the same .h kg positive increment as before. The results 
of the first three resizing steps, starting with the 270.1 m/s EAS configura- 
tion, are given in Table A-8. Comparison with Table A-5 indicates that the 
first three steps of the gradient projection search are less effective than 
the initial step of the procedure of Section A. 3. 3- The evaluation of the 
gradient projection search was terminated at this point. 






Design Variable Mass, kg 


Flutter Speed 

Step 

1 

2 

3 

4 

5 

6 

7 

8 

Total 

m/s EAS 

0 

113.5 

105.2 

90.4 

65.5 

63.6 

31.3 

21.8 

15.3 

506.6 

270.1 

1 

105.4 

97.4 

83.3 

61.0 

62.3 

34.5 

25.2 

10.3 

479.4 

270.0 

2 

98.1 

90.6 

77.1 

57-4 

61.9 

37.8 

28.2 

7.2 

458.3 

270.0 

3 

91.2 

84.1 

71.3 

54.2 

61.8 

40.7 

30.8 

5.4 

439.5 

270.0 


Table A-8: Mass Reductions with Gradient Projection Search 
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A . 5 INTERIOR PENALTY FUNCTION METHOD 


The interior penalty function method, described in Reference l 6 and 
discussed in Section 6 . 4 , employs a series of unconstrained minimizations 
of a modified objective function in order to minimize the objective function 


of interest (usually total mass). The modified objective function, 0 (m^), 
is formed by adding penalty terms, reflecting the contraint. equations, to the 
objective fiinction W(m^) (equation (A.IO)). The minimization of the modified 
objective 


0(mj^) 


W(m^) 


V r 


(A. 10 ) 


function is carried out for repeated reductions of the penalty weighting fac- 
tor, r, until the minimum of the modified objective function approximates 
the minimum of the objective function. In the present case, the penalty 
terms represented minimum size (mass) constraints for each of the eight design 
variables, in addition to the flutter speed constraint. Starting with the 
280.4 m/s EAS (545 KEAS) configuration, the initial penalty weighting factors 
were chosen to produce a total of the penalty terms approximately equal to 
the objective function. These penalty weighting factors were reduced in five 
steps to a value which resulted in penalty terms approximately equal to one 
percent of the minimim value of the objective function (total mass). The 
move direction was generated using approximate second derivatives with Newton's 
method (Section 6.4 and Reference I6) . The results of the five resizing steps 
are given in Table A- 9 • 


Step 

— 

1 

2 

3 

Desigr 

4 

Varia 

5 

ble Me 
6 

iss, kg 

7 

r 

> 

8 

Total 

Modified 

Objective 

Function 

kg 

0 

i 4 i .3 

131.0 

112.5 

81.6 

79.2 

38.9 

27.1 

19.1 

630.9 

1305.0 

1 

73.5 

75.3 

76.3 

76.9 

88.9 

52.9 

36.8 

23.4 

504.1 

1178.4 

2 

31.8 

32.8 

34.1 

47.7 

82.1 

71.0 

49.4 

14.2 

363.1 

492.6 

3 

i 4.5 

15.6 

17.7 

37.4 

81.4 

70.5 

48.7 

9.4 

295.1 

329.3 

4 

2.9 

3.5 

5.4 

33.5 

83.8 

73.4 

52.4 

8.5 

263.4 

274.4 

5 

2.0 

2.1 

1.8 

30.9 

84.1 

73.5 

52.1 

8.5 

255.1 

258.1 


Table k- 9 \ Mass Reductions with Penalty Function Procedure 
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It is recognized that additional resizing steps in which the penalty 
weighting factor is further reduced would lead to a lower total mass . How- 
ever, it has been observed (Section 6 .k. 2 ) that each step in the penalty 
function procedure requires approximately the same number of numerical evalua- 
tions as three steps in arbitrary step size procedures, so that the number of 
steps in this numerical evaluation of the penalty function method is suffi- 
cient to provide the comparison with the arbitrary step size procedures 
(Tables A-k and A- 5 ) . 


A . 6 METHOD OF FEASIBLE DIRECTIONS 

The method of feasible directions (Reference IT and Section 6 . 5 ) gen- 
erates a series of resizing move vectors, each of which is both feasible 
(does not violate active constraints) and usable (reduces total mass). Each 
resizing direction is followed until a new constraint is violated, an active 
constraint is re-encountered or the total mass is minimized. The resizing 
direction is found using the Simplex procedure to determine an optimum move 
direction. 

In the present case, the 270.1 m/s EAS configuration was taken as the 
starting point, and minimiom size (mass) constraints were imposed on the 
eight design variables in addition to the flutter speed constraint. The 
results of the first eight resizing steps are presented in Table A- 10 . 






Design Variable Mass, kg 


Flutter Speed 

Step 

1 

2 

3 

4 

5 

6 

7 

8 

Total 

m/s EAS 

0 

113.5 

105.2 

90.4 

65.5 

63.6 

31.3 

21.8 

15.3 

506.6 

270.1 

1 

98.6 

90.3 

75.5 

50.7 

69.3 

46.1 

36.7 

0.5 

467.6 

273.0 

2 

5 h .9 

46.6 

31.8 

7.0 

28.1 

89.8 

80.4 

44.1 

382.7 

270.1 

3 

23.6 

15.3 

0.5 

38.3 

59.4 

110.6 

49.0 

12.8 

309.2 

273.0 

4 

11.2 

2.9 

0.5 

47.7 

71.1 

98.0 

61.4 

0.5 

293.8 

273.7 

5 

0.5 

2.9 

0.5 

36.9 

82.5 

87.2 

61.4 

11.2 

283.0 

274.0 

5 A 

0.5 

2.9 

0.5 

30.8 

76.5 

81.2 

55.3 

5.3 

253.0 

270.1 

6 

0.5 

0.5 

0.5 

31.9 

79.0 

78.7 

52.8 

7.8 

251.6 

270.1 

6 a 

0.5 

0.5 

0.5 

31.8 

78.9 

78.5 

52.7 

7.6 

250.8 

270.1 


Table A- 10 : Mass Reductions with Feasible Directions Procedure 
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It should be noted that steps 5 A and 6 a are so designated because the 
flutter speed constraint was not active for those steps , and therefore no 
flutter speed derivatives were required. The resulting move direction was 
then one of "steepest descent," determined by the mass gradient. These two 
steps required somewhat less computing resources than did the remainder of 
the steps, and it therefore seemed reasonable not to identify them as full 
steps. It should also be noted that 0.5 kg (one pound) was established as 
the minimum mass (size) allowed for any design variable, rather than zero 
as in other procedures. This was later determined to be unnecessary as is 
discussed in Section 6 . 5 . 2 . Had this not been the case, the total mass in 
each of the last several steps would have been slightly less . 


A . 7 M OPTIMIZATION METHOD USING INCREMENTED FLUTTER ANALYSIS 

Reference l 8 describes a resizing procedure developed at Lockheed and 
used primarily in an interactive mode employing graphics displays. This 
procedure is discussed in Section 6 . 6 . 

Although this procedure is presently in an incomplete state of develop- 
ment, it is of Interest in that it employs a unique method of detennining 
resizing step size. This is done by a direct minimization of the objective 
function (weight) using a direction determined by velocity derivatives and 
the constraints. This minimization results from, and takes into account, 
the nonlinear relationship between design variable increments and flutter 
speed. A new direction is generated after each minimization, and the process 
is repeated until an acceptable approximation of the minimum total mass is 
obtained. 

The results of 4 resizing steps are presented in Table A- 11 . 






Design Variable Lfess, kg 


Flutter Speed 

step 

1 

2 

3 

4 

5 

6 

7 

8 

Total 

m/s EAS 

0 

113.5 

105.2 

90.4 

65.5 

63.6 

31.3 

21.8 

15.3 

506.6 

270.1 

1 

0.0 

0.0 

0.0 

70.1 

71.6 

44.0 

63.5 

19.5 

268.6 

270.1 

2 

0.0 

0.0 

0.0 

40.3 

84.8 

78.2 

49.4 

0.0 

251.8 

270.1 

3 

0.0 

0.0 

0.0 

36.0 

82.7 

73.2 

47.7 

10.8 

250.3 

270.1 

4 

0.0 

0.0 

0.0 

32.7 

84.9 

72.7 

52.9 

6.7 

249.8 

270.1 


Table A-ll: Mass Reductions with the Method Using 

Incremented Flutter Analysis 
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